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ABSTRACT 

The relation between galaxy stellar mass and gas-phase metallicity is a sensitive diagnostic of the main 
processes that drive galaxy evolution, namely cosmological gas inflow, metal production in stars, and gas 
outflow via galactic winds. We employed the direct method to measure the metallicities of ^200,000 star- 
forming galaxies from the SDSS that were stacked in bins of (1) stellar mass and (2) both stellar mass and star 
formation rate (SFR) to significantly enhance the signal-to-noise ratio of the weak [O III] A4363 and [O II] 
AA7320, 7330 auroral lines required to apply the direct method. These metallicity measurements span three 
decades in stellar mass from log(Mi/M Q ) = 7.4-10.5, which allows the direct method mass-metallicity relation 
to simultaneously capture the high-mass turnover and extend a full decade lower in mass than previous studies 
that employed more uncertain strong line methods. The direct method mass-metallicity relation rises steeply 
at low mass (O/H oc M* I/2 ) until it turns over at log(M*/M Q ) = 8.9 and asymptotes to 12 + log(0/H) = 8.8 
at high mass. The direct method mass-metallicity relation has a steeper slope, a lower turnover mass, and a 
factor of two to three greater dependence on SFR than strong line mass-metallicity relations. Furthermore, the 
SFR-dependence appears monotonic with stellar mass, unlike strong line mass-metallicity relations. We also 
measure the N/O abundance ratio, an important tracer of star formation history, and find the clear signature of 
primary and secondary nitrogen enrichment. N/O correlates tightly with oxygen abundance, and even more so 
with stellar mass. 

Subject headings: Galaxies: general — Galaxies: abundances — Galaxies: ISM — Galaxies: evolution — 
Galaxies: stellar content — ISM: abundances 



1. INTRODUCTION 

Galaxy metallicities are one of the fundamental observa- 
tional quantities that provide information about their evolu- 
tion. The metal content of a galaxy is governed by a complex 
interplay between cosmological gas inflow, metal production 
by stars, and gas outflow via galactic winds. Inflows dilute 
the metallicity of a galaxy in the short term but provide the 
raw fuel for star formation on longer timescales. This gas 
turns into stars, which convert hydrogen and helium into heav- 
ier elements. The newly formed massive stars inject energy 
and momentum into the gas, driving large-scale outflows that 
transport gas and metals out of the galaxy. The ejected metals 
can escape the gravitational potential well of the galaxy to en- 
rich the intergalactic medium or reaccrete onto the galaxy and 
enrich the inflowing gas. This cycling of baryons in and out 
of galaxies directly impacts the stellar mass (M*), metallicity 
(Z), and star formation rate (SFR) of the galaxies. Thus, the 
galaxy stellar mass-metallicity relation (MZR) and the stellar 
mass-metallicity-SFR relation serve as crucial observational 
constraints for galaxy evolution models that attempt to under- 
stand the build up of galaxies across cosmic time. Here we 
present new measurements of the MZR and the M^-Z-SFR 
relation that span three orders of magnitude in stellar mass 
with metallicities measured with the direct method. 

The first indication of a correlatio n betw een mass and 
metallicity came when iLequeux et al.1 (119791) demonstrated 
the existence of a relation between total mass and metallicity 
for irregular and blue compact galaxies. Subsequent studies 
showed that metallicity als o correlates with o ther galaxy prop- 
erties, such as luminosity (Rubin et al. 1984) and rotation ve- 
locity (IZaritsky et al.lll994t lGarnettll2002l) . The advent of re- 
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liable stellar population synthesis models dBruzual & CharloU 
120031) enabled more accurate ste llar mass measu r ement s from 
spectral energy distributions. iTremonti et al.l (12004. here- 
after T04) showed the existence of a tight correlation between 
galaxy stellar mass and metallicity among ~53,000 galaxies 
from t he Sloan Digital Sky Survey (SDSS; lYork et al.ll2000l) 
DR2 ( Abazaiian et al. 2003) base d on th e tell ar m ass mea- 
surements from lKauffmann et al.l (l2003al) . The HOI MZR in- 
creases as roughly O/H ex M* 1 / 3 fromM* = 10 8 5 -10 10 5 M© 
and then flattens above M* ~10 105 M Q . They found that 
the scatter in the MZR was smaller than the scatter in the 
luminosity-metallicity relation an d concluded tha t the MZR 
was more physically motivated. iLee et all (120061) extended 
the MZR down another ^2.5 dex in stellar mass with a sam- 
ple of local dwa r f irreg ular galaxies. The scatter and s lope 
of the ILee et al.l (120061) MZ R are c onsistent with th e IT04I 
MZR (cf., Zahi d et al.ll2012ah . but the lLee et alJ J2006I) MZR 
is offset to lower metallicitie s by 0.2-0.3 de x. This off- 
set is likely because T04 and Lee et al.l (120061) use different 
metho ds to estimate metallicity. Later work by lEllison et al.1 
(2008) discovered that galaxies with high SFRs (and larger 
half-light radii) are systematically offset to lower metallici- 
ties than more weakly star-form ing galaxies at the same stel- 
lar mass. iMannucci et al.l (120101) and lLara-Lopez et al.l (120101) 
studied this effect in a systematic fashion and demonstrated 
that the scatter in the MZR is r educed further by accounting 
for SFR. Mann ucci et al.l (120101) introduced the concept of the 
fundamental metallicity relation (FMR) by parametrizing the 
second-order dependence of the MZR on SFR with a new ab- 
scissa, 

Ma =log(^)-alog(SFR), (1) 

where the coefficient a is chosen to minimize the scatter in 
the relation. We will refer to this particular parametriza- 
tion as the FMR but the general relation as the M+-Z- 
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SFR relation. I nteres tingly, iMannucci et al.l ( 120101) and 
lLara-Lopez et al.1 (120101) found that the M^-Z-SFR relation 
does n ot evolve with redshift up to z y 2.5 , as opposed to the 
MZR dErb et alj2006tlMaiolino et alj2"008tlZahid et al.l201 It 
Moustakas et al .l201 lh . However, this result depends on chal- 
len ging high redshi ft metallicity measurements, specifically 
the Er b et ajj (12006b sample of sta cked galaxy spectra at z ~ 
2.2 and the Maiolin o et alj Q2008) sample of nine galaxies at 
z~3.5. 

Galaxy evolution models aim to reproduce various features 
of the MZR and M+-Z-SFR relation, specifically their slope, 
shape, scatter, and evolution. The most distinguishing charac- 
teristic of the shape of the MZR is that it appears to flatten and 
become independent of mass at ~ 10 10 5 M Q . The canon- 
ical explanation is that this turnover reflects the efficiency of 
metal ejection from galaxies because the gravitational poten- 
tial wells of galaxies at and above this m ass scale are too deep 
for supernova-driven winds to escape jPekel & Silkl 1986: 
iDekel & Wooll2003t iTremonti et al. 1 120041) . In this scenario, 
the metallicity of these galaxies approaches the effective yield 
of the stellar populat i on. H o wever, recent simulatio ns by 
Qppenheimer & Pavel (120061) . iFinlator & Pavel (120081) . and 
Pave et all (201 la b) show that winds characterized by a con- 
stant velocity and constant mass-loading parameter (mass 
outflow rate divided by SFR; their cw simulations), which 
were intended to represent supernova-driven winds, result in 
a MZR that fails to qualitatively match observations. The cw 
simulations produce a MZR that is flat with a very large scat- 
ter at low mass, yet becomes steep above the blowout mass, 
which is the critical scale above which all metals are retained. 
Instead, they find that their simulations with momentum- 
driven winds dMurray et al.ll2005t IZhang & Thompsonll2012l) 
best reproduce the slope, shape, scatter, and evolution of the 
MZR because the wind velocity scales with the escape veloc- 
ity of the halo. Their model naturally produces a FMR that 
show s little evolution since z = 3, consistent with observa- 
tions (IMannucci et al.ll20l61 : [Richard et al. 2 0"Tlt ICresci et alj 
120121) . However, their FM R does not quite reach the low ob- 
served scatter reported by IMannucci et al.l (120101) . Addition- 
ally, they find that the coefficient relating and SFR that 
minimiz es the scatter in the F MR is different from the one 
found by IMannucci et alj d2010t) . While there is hardly a con- 
sensus among galaxy evolution models about how to produce 
the MZR and M^-Z-SFR relation, it is clear that additional 
observational constraints would improve the situation. So far, 
the overall normalization of the MZR and the M*-Z-SFR re- 
lation have been mostly ignored by galaxy evolution models 
due to uncertainties in the nucleosynthetic yields used by the 
models and the large (up to a factor of five) uncertainties in the 
normalization of the observed relations caused by systematic 
offsets among metallicity calibrations. If these uncertainties 
could be reduced, then the normalization could be used as an 
additional constraint on galaxy evolution models. 

The current metallicity and the metal enrichment history 
also have implications for certain types of stellar explosions. 
There is mounting evide nce that long duration gamma ray 
bursts ( Stanek et al. 2006), over-luminous type II supernovae 
(IStoll et al.l 120111). and super-Chandrasekhar type la super- 
novae (Khan et al. 201 1) preferentially occur in low metallic- 
ity environments. The progenitors of long gamma ray bursts 
and over-luminous type II supernovae are thought to be mas- 
sive stars and the nature of their explosive death could plau- 
sibly depend on their metallicity. The cause of the associ- 



ation between super-Chandrasekhar type la supernovae and 
low metallicity environments is still highly uncertain because 
the progenitors are not well known. Nevertheless, accurate 
absolute metallicities for the host galaxies of the progenitors 
of gamma ray bursts, over-luminous supernovae, and super- 
Chandrasekhar type la supernovae will help inform the mod- 
els of stellar evolution and explosions that attempt to explain 
these phenomena. 

The uncertainty in the absolute metallicity scale can be 
traced to differences between the two main methods of mea- 
suring metallicity: the direct method and strong line method. 
The direct method utilizes the flux ratio of auroral to strong 
lines to measure the electron temperature of the gas, which is 
a good proxy for metallicity because metals are the primary 
coolants of H II regions. This flux ratio is sensitive to tem- 
perature because the auroral and strong lines originate from 
the second and first excited states, respectively, and the rela- 
tive level populations depend heavily on electron temperature. 
The electron temperature is a strong function of metallicity, 
such that hotter electron temperatures correspond to lower 
metallicities. In the direct method, the electron temperature 
estimate is the critical step because the uncertainty in metal- 
licity is nearly always dominated by the uncertainty in the 
electron temperature. The strong line method uses the flux 
ratios of the strong lines, which do not directly measure the 
metallicity of the H II regions but are metallicity-sensitive and 
can be calibrated to give approximate metallicities. The direct 
method is chosen over strong line methods when the auroral 
lines can be detected, but these lines are often too weak to de- 
tect at high metallicity. The strong lines, on the other hand, are 
much more easily detected than the auroral lines, particularly 
in metal-rich objects. Consequently, the strong line method 
can be used across a wide range of metallicity and on much 
lower signal-to-noise ratio (SNR) data, so nearly all metal- 
licity studies of large galaxy samples employ the strong line 
method. Pespite the convenience of the strong line method, 
the relationship between strong line ratios and metallicity is 
complicated due to the sensitivity of the strong lines to the 
hardness of the incident stellar radiation field and the exci- 
tation and ionization states of the gas. Thus, strong line ra- 
tios must be calibrated (1) empirically with direct method 
metallicities, (2) theoretically with photoionization models, 
or (3) semi-empirically with a combination of direct method 
metallicities and theoretically calibrated metallicities. Unfor- 
tunately, the three classes of calibrations do not genetically 
produce consistent metallicities. For example, metallicities 
determined with theoretical strong line calibrations are sys- 
tematically higher than those from the direct method or em- 
pirical strong line c alibrations by up to ~0.7 dex (for a de - 
tailed discussion see Moustakas et al. 2010; Stasihska 2010). 
The various strong line methods also exhibit systematic dis- 
agreements as a function of metallicity and perform better or 
poorer in certain metallicity ranges. 

The cause of the discrepancy between direct method metal- 
licities and theoretically calibrated met a llicitie s is currently 
unknown. As recognized by Peimbert (119671) . the electron 
temperatures determined in the direct method might be over- 
estimated in the presence of temperature gradients and/or fluc- 
tuations in H II regions. Such an effect wo uld cause the di- 
rect method m etallicities to be biased low ( Stasihska! 120051 ; 
iBresohnl 1200 8). A similar result could arise if the tradition- 
ally adopted electron energy distrib ution is different fro m the 
true distribution, as suggested by Nich olls et al.l (120121) . Al- 
ternatively, the photoionization models that serve as the basis 
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for the theoretical st rong line calibratio ns, such as CLOUDY 
dFerland et al.l [1998) and mappings dSutherland & Dopital 
119931) . make simplifying assumptions in their treatment of 
H II regions that may result in overestimated metallicities, 
such as t he geometry of the ne bula or the age of the ionizing 
stars (see Moustakas et al. 2010, for a thorough discussion of 
these issues); however, no one particular assumption has been 
conclusively identified to be the root cause of the metallicity 
discrepancy. 

In this work, we address the uncertainty in the absolute 
metallicity scale by using the direct method on a large sample 
of galaxies that span a wide range of metallicity. The uniform 
application of the direct method also provides more consis- 
tent metallicity estimates over a broad range in stellar mass. 
While the auroral lines used in the direct method are unde- 
tected in most galaxies, we have stacked the spectra of many 
galaxies (typically hundreds to thousands) to significantly en- 
hance the SNR of these lines. In Section [2] we describe 
the sample selection, stacking procedure, and stellar contin- 
uum subtraction. Section [3] describes the direct method and 
strong line metallicity calibrations that we use. In Section|4j 
we demonstrate that mean galaxy properties can be recovered 
from stacked spectra. We show the e lectron temperature re- 
lations for the stacks in Section |3~T1 and argue that r e [0 II] 
is a bette r tracer of oxygen abundance than T e [0 III] in Sec- 
tion [52] Section [5] shows the main results of this study: the 
MZR and M+-Z-SFR relation with the direct method. In Sec- 
tion [6] we present the direct method N/O relative abundance 
as a function of O/H and stellar mass. Section [7] details the 
major uncertainties in metallicity measurements and the im- 
plications for the physical processes that govern the MZR and 
M+-Z-SFR relation. Finally, we present a summary of our re- 
sults in Section |8] For the purpose of discussing metallicities 
relative to the solar value, we ado pt the solar oxygen abun- 
dance of 12 + log(0/H) = 8.86 from Delahave & Pinsonneaulll 
(2006). Throughout this work, stellar masses and SFRs are 
in units of M Q and M Q yr" 1 , respectively. We assume a 
standard ACDM cosmology with O m = 0.3, J7a = 0.7, and 
// = 70kms" 1 Mpc" 1 . 

2. METHOD 
2.1. Sample Selection 

The observations for our g alaxy sample come f rom the 
SDSS Data Release 7 (DR7; Abazaiia netai] 120091). a s ur- 
vey that includes -930,000 galaxies dStrauss et al.l 120021) in 
an area of 8423 square degrees. The parent sample for this 
study comes from the MPA-JHU catalog o f 818,333 unique 
galaxie s which have derived stellar masses ([Kauffman n et al] 
I2003al) . SFRs (Brinchm ann et al.l l2004t ISalim et all 12007b . 
and metallicities (IT04I) . We chose only galaxies with reliable 
redshifts (er- < 0.001) in the range 0.027 < z < 0.25 to ensure 
that the [O II] A3727 line and the [O II] AA7320, 7330 lines 
fall within the wavelength range of the SDSS spectrograph 
(3800-9200 A). 

We discard galaxies classified as AGN because AGN emis- 
sion line ratios may produce spurious metal licity measure- 
ments. We adopt the Kauffm ann et alj d2003bl) criteria (their 
Equation 1) to differentiate between star-forming galaxies and 
AGN, which employs the emission line ratios that define the 
Baldwin, Phillips, and Terlevich (1981) (BPT) diagram: 

log([0 in] A5007/H/3) > 

2 Available at http://www.mpa-garching.mpg.de/SDSS/DR7/ 



0.61 [log([N II] A6583 /Ha) - 0.05] -1 + 1.3. (2) 

We follow thellOlSNR thresholds for emission lines. Specif- 
ically, we restrict our sample to galaxies with H/3, Ha, and 
[N II] A6583 detected at >5a. Further, we apply the AGN- 
star-forming galaxy cut (Equation [2j to galaxies with >3a 
detections of [O III] A5007. We also select galaxies with 
[O in] A5007 < 3ct but log([N n] A6583/Ha) < -0.4 as 
star-forming to include high metallicity galaxies with weak 
[O in] A5007. 

At the lowest stellar masses (logfM*] < 8.6), this ini- 
tial sample is significantly contaminated by spurious galax- 
ies, which are actually the outskirts of more massive galax- 
ies and were targeted due to poor photometric deblending. 
We remove galaxies whose photometric flags include DE- 

BLEND_NOPEAK or DEBLENDED_AT_EDGE. We also visu- 
ally inspected all galaxies with log(M*) < 8.6 and discarded 
any that suffered from obvious errors in the stellar mass de- 
termination (again, likely as a result of off-center targeting of 
a much more massive galaxy). 

After all of our cuts, the total number of galaxies in our 
sample is 208,529 and the median redshift is z = 0.078. At 
this redshift, the 3" diameter SDSS aperture will capture light 
from the inner 2.21 kpc of a galaxy. Since t he central r egions 
of galaxies will tend to be more metal-rich ( Searlel fl 97 lh . the 
metallicities measured from these observations will likely be 
biased high due to the aperture size relative to angular extent 
of the galaxies. However, we expect this bias i s small for most 
galax i es (for a more deta iled discussion see iTremonti et al.l 
l2004t iKewley et alJl2005l) . In particular, the galaxies with 
very low stellar masses and metallicities that define the low 
mass end of the MZR ten d to be compact and have ho moge- 
neous metallicities (e.g., iKobulnickv & Skillmanlll997T) . al- 
th ough many of these are excluded by the criteria proposed 
by lKewley etaP J2005I) . 

2.2. Stacking Procedure 

The primary motivation for this investigation is to mea- 
sure the metallicity of galaxies with the direct method. The 
main challenge is that the weak [O III] A4363 and [O II] 
AA7320, 7330 auroral lines are undetected in most of the in- 
dividual spectra. To improve the SNR of the spectra, we 
stacked galaxies that are expected to have similar metallici- 
ties and hence line ratios. Given the tightness of the MZR 
and M^-Z-SFR relation, it is reasonable to expect that galax- 
ies at a given stellar mass, or simultaneously a given stellar 
mass and SFR, will have approximately the same metallicity. 
Thus, we have created two sets of galaxy stacks: (1) galax- 
ies binned in 0.1 dex in M* from log(MJM Q ) = 7.0 to 11.0 
(hereafter stacks) and (2) galaxies binned in 0. 1 dex in 
M* from log(M*/M Q ) = 7.0 to 1 1.0 and 0.5 dex in SFR from 
log(SFR/[M Q yr" 1 ]) = -2.0 to 2 (hereafter M^-SFR sta cks). 
We adopt the total stellar mass ([Kau ffmann et al. 2003 ah and 
the total SFR dBrinchmann et al]|2004t ISalim et al.ll2007l) val- 
ues from the MPA-JHU catalog, as opposed to these quanti- 
ties calculated only for the light within the fiber. For conve- 
nience, we will refer to the stacks by the type of stack with 
a subscript and a superscript to denote the upper and lower 
bounds of log(M^) or log(SFR) (e.g., M*|;| is the M* stack 
with log[M*/M Q ] = 8.7-8.8, and SFR[J;^ corresponds to the 
M*-SFR stacks with log[SFR/M Q yr -1 ] = 0.0-0.5). Figure 
[T] shows the number of galaxies in each M^-SFR stack (each 
box represents a stack) with a measured metallicity (indicated 
by the color coding). 
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Figure 1. Number of galaxies and direct method metallicity as a function of M* and SFR. The squares represent each M*-SFR stack, the number of galaxies is 
indicated by the white text, and the color scale corresponds to the metallicity. For reference, the Tremonti et al. (2004) MZR covers log(M*) = 8.5-1 1.5, and the 
IMannucci e t al] <2O10l) FMR spans log(M*) = 9.1-1 1.35 and log(SFR) = -1.45-^0.80. 
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Figure 2. Sample spectra from the log(M*) = 8.7-8.8 (AT ga i = 884) stack. From left to right, the three columns show the [O III] A4363, [N II] A5755, and [O II] 
AA7320, 7330 auroral lines. From top to bottom, the four rows correspond to the reduced spectrum of a single galaxy, the spectrum of the stack, the spectrum 
of the stack after the removal of the stellar continuum (fit from 3700-7360 A), and the spectrum of the stack after the removal of the stellar continuum (fit to a 
200 A window near the emission line of interest). The continuum rms of each spectrum near the relevant emission line is given in the inset of each panel. 



We stacked galaxy spectr a that have been proce ssed with 
the SDSS reduction pipeline (Stou ghton et al .120021) . First, we 
corre cted for Milky Way r eddening with the extinction values 
from|Schlegel et al. ( 1998). Then, the individual galaxy spec- 
tra were shifted to the rest frame with the redshifts from the 
MPA/JHU catalog. Next, we linearly interpolated the spectra 
onto a universal grid (3700-7360 A; AA = 1 A) in linear-A 
space. This interpolation scheme conserves flux in part be- 
cause the wavelength spacing of the grid is narrower than the 
width of bright emission lines. The spectra were then normal- 
ized to the mean flux from 4400^-450 A. Finally, the spectra 
were co-added (i.e., we took the mean flux in each wavelength 
bin) to form the stacked spectra (see Section |4] for compar- 
isons between the electron temperatures and metallicities of 
stacks and individual galaxies). 

Figure [2| shows the SNR increase of the [O III] A4363 
(left column), [N II] A5755 (middle column), and [O II] 
AA7320, 7330 (right column) lines as the spectra are pro- 
cessed from a typical single galaxy spectrum (top row) to 
the stacked spectrum (second row) to the st ellar continuum 
subtracted spectrum (third row; see Section 12.3) or the nar- 
row wavelength window stellar continuum subtracted spec- 



trum (bottom row; see Section 12.3b . The spectra in the top 
row are from a typical galaxy in the log(M^) = 8.7-8.8 bin; 
the bottom three rows show the stacked spectra from the same 
bin. In each panel, we report the continuum root mean square 
(rms). The decrease in the continuum noise when comparing 
the spectra in the top row to the second row of Figure |2]is 
dramatic. Further significant noise reduction can be achieved 
by removing the stellar continuum (shown i n the bottom two 
rows of Figure|2), as we describe in Section |2~31 

2.3. Stellar Continuum Subtraction 

Stacking the spectra increases the SNR, but it is im- 
portant to fit and subtract the stellar continuum to detect 
and accurately measure the flux of these lines, especially 
[O III] A4363 due to its proximity to the H7 stellar absorp- 
tion feature. We subtracted the stellar continuum with syn- 
thetic template galaxy spectra created with the S TARLIGHT 
stel lar synthesis code (ICid Fernandes et al.l 120051) . adopted 
the iCardelli et al. I (119891) extinction law, and masked out 
the locations of the emission lines. The synthetic spectra 
were created f rom a library of 300 empirica l MILES spec- 
tral templates (Sanchez-Blazquez et al. 20061: ICenarro et all 
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Table 2 

Line Fluxes 



Column Format Description 



Table 1 

Wavelength Fit and Mask Ranges of 
Measured Lines 



Line 



(1) 



Fit Range 
[A] 

(2) 



Mask Range 
[A] 

(3) 



[O II] A3727 
[Ne III] A3868 
[S II] A4069 

H7 A4340 
[O III] A4363 
He II A4686 
[Ar IV] A4740 

A4861 
[O in] A4959 
[O III] A5007 
[N II] A5755 
[S III] A6312 
[N II] A6548 

Ha A6563 
[N II] A6583 
[S II] A6716 
[S II] A6731 
[Ar III] A7135 
[O II] A7320 
[O II] A7330 



3700-4300 
3800-4100 
3950-4150 
4250-4450 
4250-4450 
4600-4800 
3700-7360 
3700-7360 
3700-7360 
3700-7360 
5650-5850 
6100-6500 
3700-7360 
3700-7360 
3700-7360 
3700-7360 
3700-7360 
7035-7235 
7160-7360 
7160-7360 



3710-3744 
3863-3873 

4336-4344 
4360^1366 
4680^1692 

4857^1870 
4954^964 
5001-5013 
5753-5757 
6265-6322 
6528-6608 
6528-6608 
6528-6608 
6696-6752 
6696-6752 
7130-7140 
7318-7322 
7328-7332 



Note. — Column (1): Emission lines. 
Column (2): The wavelength range of the 
stellar continuum fit. Column (3) The 
wavelength range of the stellar contin- 
uum fit that was masked out. 
l2007t iVazdekis etafll20Tot iFalcon-Barroso et al.ll20TTl data 
as obtained from the MILES website0). The MILES tem- 
plates provided an excellent fit to the stellar continuum (see 
bottom two rows of Figure |2). We note the MILES tem- 
pla tes yielded better fits to t he very high SNR spectra than 
the Bruzual & Chariot (2003) spectral templates, based on the 
STELIB (Le Borgne et al. 2003) library. 

We performed stellar template fits to the entire spectral 
range, select subregions centered on weak lines of interest, 
and subregions around the strong lines blueward of 4000 A. 
The latter are situated among a forest of stellar absorption 
lines. The line fluxes of the strong emission lines redward of 
4000 A (H/3, [O in] AA4959, 5007, Ha, [N 11] AA6548, 6583, 
and [S II] AA6716, 6731) were measured from the spectrum 
where the stellar continuum was fit over the full wavelength 
range of our stacked spectra ( A = 3600-7360 A; see third row 
of Figure The stellar continuum subtraction near weak 
emission lines ([S II] A4069, [O III] A4363, He II A4686, 
[N 11] A5755, [S in] A6312, [Ar iv] A4740, and [O 11] 
AA7320, 7330) and blue strong emission lines ([O II] A3727 
and [Ne III] A3868) was improved if the stellar continuum 
fit was restricted to limited wavelength ranges within a few 
100 A of the line of interest (compare the third and bottom 
rows of Figure|2]i. For the weak lines and blue strong lines, we 
measured the line fluxes from the stellar continuum subtracted 
spectra within these narrow wavelength windows (details are 
listed in Table[TJ. In order to compare the line fluxes across re- 
gions with different stellar continuum subtraction (e.g., from 
portions of the spectrum that were fit with smaller wavelength 
ranges), we denormalized the spectra after the STARLIGHT fit. 

2.4. Automated Line Flux Measurements 
3 http://miles.iac.es/ 
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F4.1 


Lower stellar mass limit of the stack 
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F4.1 


Upper stellar mass limit of the stack 
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F4.1 


Lower SFR limit of the stack 
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F4.1 


Upper SFR limit of the stack 
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15 


Number of galaxies in the stack 
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F6.3 


Median stellar mass of the stack 
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F6.3 


Median SFR of the stack 
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F6.2 


[O II] A3727 line flux 
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F5.2 


Error on [O II] A3727 line flux 


10 


F5.2 


[Ne III] A3868 line flux 


11 


F4.2 


Error on [Ne III] A3868 line flux 


12 


F4.2 


[S II] A4069 line flux 


13 


F4.2 


Error on [S II] A4069 line flux 


14 


F6.2 


H7 line flux 


15 


F5.2 


Error on H7 line flux 


16 


F5.2 


[O III] A4363 line flux 


17 


F4.2 


Error on [O III] A4363 line flux 


18 


F4.2 


He II A4686 line flux 


19 


F4.2 


Error on He II A4686 line flux 


20 


F4.2 


[Ar IV] A4740 line flux 


21 


F4.2 


Error on [Ar IV] A4740 line flux 


22 


F6.2 


[O III] A4959 line flux 


23 


F4.2 


Error on [O III] A4959 line flux 


24 


F6.2 


[O III] A5007 line flux 


25 


F5.2 


Error on [O III] A5007 line flux 


26 


F4.2 


[N II] A5755 line flux 


27 


F4.2 


Error on [N II] A5755 line flux 


28 


F4.2 


[S III] A63121ine flux 


29 


F4.2 


Error on [S III] A6312 line flux 


30 


F5.2 


[N II] A6548 line flux 


31 


F4.2 


Error on [N II] A6548 line flux 


32 


F6.2 


Ho line flux 


33 


F5.2 


Error on Ha line flux 


34 


F6.2 


[N II] A6583 line flux 


35 


F4.2 


Error on [N II] A6583 line flux 


36 


F6.2 


[S II] A6716 line flux 


37 


F4.2 


Error on [S II] A6716 line flux 


38 


F5.2 


[S II] A6731 line flux 


39 


F4.2 


Error on [S II] A6731 line flux 


40 


F4.2 


[Ar III] A7135 line flux 


41 


F4.2 


Error on [Ar III] A7135 line flux 


42 


F4.2 


[O II] A7320 line flux 


43 


F4.2 


Error on [O II] A7320 line flux 


44 


F4.2 


[O II] A7330 line flux 


45 


F4.2 


Error on [O II] A7330 line flux 



Note. — This table is published in its entirety in the 
electronic edition of the journal. The column names 
are shown here for guidance regarding its form and 
content. 

We used the specfit task dKrisslll994h in the IRAF/STSDAS 
package to automatically fit emission lines with a x 2 mini- 
mization algorithm. We simultaneously fit a flat continuum 
and Gaussian line profiles for the emission lines, even if 
lines were blended. For doublets, we fixed the width of the 
weaker line by pinning its velocity width to the stronger line 
([O 11] A3726 to [O 11] A3729, [O m] A4959 to [O m] A5007, 
[N 11] A6548 to [N 11] A6583, [S 11] A6731 to [S 11] A6716, 
and [O 11] A7330 to [O 11] A7320). We also included the 
continuum rms of the spectrum as an input to the fitting pro- 
cedure. After experimenting with several different x 2 min- 
imization algorithms implemented within specfit, we chose 
the simplex algorithm because of its consistent convergence, 
particularly for weak lines. Line fluxes measured by spec- 
fit generally agreed well with line fluxes measured interac- 
tively with the OSU LINER package. The uncertainty in the 
line flux is based on the \ 2 fit returned from specfit. Fi- 
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nally, all line fluxes were correcte d for r eddening with the 
extinction law from iCardelli et al.l (119891) and the assump- 
tion that the intrinsic ratio of the Balmer lines is set by case 
B recombination (Ha/H/3 = 2.86 for T e = 10,000 K). We 
adopted a fixed Ha/H/3 ratio, even though it is a weak func- 
tion of electron temperature. For the log(M + /M Q ) = 10.0- 
10.1 stack (T e [0 II] = 7200 K), whose oxygen abundance 
is dominated by + (i.e., a stack where the potential effect 
would be maximal due to the long wavelength baseline be- 
tween [O II] A3727 and [O n] AA7320, 7330), this effect 
would decrease log(0 + /H + ) by ^0.07 dex. The line fluxes 
are presented in an online table whose columns are described 
in Table [2] 

We disregarded lines that were poorly fit (negative flux, un- 
certainty in central wavelength > 1 A, had uncertainty in the 
velocity width of >100 km/s, or had low SNR [<5cr]). Fur- 
ther care was taken to ensure the robustness of [O III] A4363 
flux measurements. As increased to moderate values 
(log[M*] > 9.0), an unidentified emission feature at 4359 A 
became blended with the [O III] A4363 line, which limited 
the SNR of the line flux measurement independent of the con- 
tinuum rms. We are unsure of the origin of this feature, but 
it could be caused by an over-subtraction in the stellar con- 
tinuum fit. We simultaneously fit the 4359 A feature and 
[O III] A4363 and pinned the velocity width of both lines 
to H7. If 4359A > 0.5 [O III] A4363, then we determined 
that [O III] A4363 could not be robustly fit. If [O III] A4363 
could be well fit, we refit it with a single Gaussian whose 
velocity width was pinned to H7. The line flux measure- 
ments from the single Gaussian fitting agreed better with 
interactive line flux measurements than the deblended line 
flux measurements. The remaining weak lines are in regions 
without strong stellar absorption features. Often, the [O II] 
AA7320, 7330 lines could be detected in the stacked spec- 
tra without the stellar continuum fit (see Figure |2f). The 
[N II] A5755 and [S II] A4069 auroral lines were usually too 
weak to be detected without stellar continuum subtraction. 

Optical recombination lines, such as C II A4267 and 
O II A4649, are also sensitive to metallicity. Unlike auro- 
ral lines, they are almost independent of temperature, so they 
could provide a useful check on the direct method metallic - 
ities. Unfortunately, optical recombination lines tend to be 
very weak (e.g., the median O II A4649/[ Q Hi] A4363 ratio of 
five extragalactic H II regions studied by Esteban et al. 2009 
was 0.08), and we did not detect them in the stacked spectra. 

3. ELECTRON TEMPERATURE AND DIRECT ABUNDANCE 
DETERMINATION 

3.1. Electron Temperatures 

Different ionic species probe the temperature of different 
ioniza tion zones of H II regions (e.g., Stasinska 1982; Garnetj 
119921) . In the two -zone model, the high ionization zone is 
traced by [O III], and t he low ionization zon e is traced by 
[O II], [N II], and [S III. Cam pbell etall (119861) used the pho- 
toionization models of lStasinskal {T982) to derive a linear re- 
lation between the temperatures in these zones, 

r e [o 11] = r e [N 11] = r e [s n] = o.7r e [o m] + 3000, (3) 

where T e is in units of K. Subsequen tly, we will refe r to 
this relation as th e T2-T3 relation (see Pa gel et al.lfT992l and 
llzotov et"ai]l2006l for alternative formulations of the T2-T3 re- 
lation). This relation is especially useful to infer the abun- 
dance of unseen ionization states, a critical step in measuring 



the total oxygen abundance. While convenient, this theoret- 
ical relation may be one of the biggest uncertainties in the 
direct method because it is not definitively constrained by ob- 
servations due to the larg e random errors in the flux of [O II] 
AA7320, 7330 (e.g., see iKennicutt et al]|2003t iPilvugin et ail 
2006). The high SNR of our stacked spectra enables us to 
measure the electron temperature of both the high and low 
ionization zones for many of our stacks. 

We measured the electron temperature of [O III], 
|"0 II], [N II], and [S 11] with the nebular.temden routine 
(IShaw & Dufourlll995l) in ir af/stsdas, which is bas ed on 
the five level atom program of De Robertis et al. ( 1987). This 
routine determines the electron temperature from the flux ratio 
of the auroral to strong emission line(s) for an assumed elec- 
tron density. The diversity of these temperature diagnostics 
are valuable cross-checks and provide an independent check 
on the applicability of the relation; however, for mea- 

suring oxygen abundances, we only use r e [0 III] and TJO II]. 
The electron density (n e ) can be measured f rom the density 
sensit ive [S 11] AA6716, 6731 doublet (cf.. ICai & Pradhanl 
[1993b . For 6/45 of the M* stacks and 65/228 of the M*-SFR 
stacks, [S 11] A6716 / [S n] A6731 was abov e the theoreti- 
cal maximum ratio of 1 .43 dOsterbrocklll989l) . which firmly 
places these galaxies in the low densit y regime, and we as- 
sume n e = 100 cm" 3 for our analysis. lYin et al.l d2007) found 
similar inconsistencies between the theoretical maximum and 
measured flux ratios for individual galaxies, which suggests 
that there might be a real discrepancy between the maximum 
observed and theoretical values of [S II] A6716 / [S II] A6731. 

We calculated the electron temperature and density uncer- 
tainties by propagating the line flux uncertainties with Monte 
Carlo simulations. For the simulations, we generated 1,000 
realizations of the line fluxes (Gaussian distributed accord- 
ing to the la uncertainty) and processed these realizations 
through nebular.temden. The electron temperatures of the 
stacks are given in Table[3](full version available online). 

In Figure [3] we plot the electron temperatures of [O II], 
[N II], and [S II] against the [O III] electron temperature for 
the M* stacks (left column; open circles) and the M*-SFR 
stacks (right column; circles color-coded by SFR). For com- 
parison, we show the To-Tt, relation (Equation[3]l as the black 
line in each panel. In all three T e -T e plots, the M t stacks 
form a tight locus that falls within the distribution of M^-SFR 
stacks. The M^-SFR stacks show a large dispersion in T e [O II] 
at fixed r e [0 III] that is not present in the M t stacks. Most of 
this scatter is due to stacks with SFRj q, which approach and 
exceed the T e [0 ll]-r e [0 III] relation. On the other hand, the 
M*-SFR stacks show little scatter in the T e [N ll]-r e [0 III] and 
r e [S ll]-r e [0 III] plots, and they track the M* stacks in these 
plots. 

The vast majority of the stacks in Figure[3]fall below the T%— 
73 relation, independent of the type of stacks (M* or M+-SFR) 
or the tracer ion ([O II], [N II], or [S II]). The multiple tem- 
perature indicators show that the T2-T3 relation overpredicts 
the temperature in the low ionization zone (or underpredicts 
the temperature in the high ionization zone). If we assume 
that r e [0 III] is accurate (i.e., the temperature in the low ion- 
ization zone is overestimated by the T2-T3 relation), then the 
median offsets from the T2-T3 relation for the M t stacks and 
the M*-SFR stacks, respectively, are 

• r e [0 II]: -2000 K and -1300 K, 

• r e [N II]: -1200 K and -1400 K, 
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Table 3 

Electron Temperatures, Metallicity, and N/O Abundance 



log(M.) 


log(SFR) 


Ngal 


7i[Om] 


T e [0 II] 


7i [Nil] 


Te[S II] 


12 + log(0/H) 


log(N/0) 


[Mq] 


[M yr" 1 ] 




[K] 


[K] 


[K] 


[K] 


[dex] 


[dex] 


(1) 


(2) 


(3) (4) 


(5) 


(6) 


(V) 


(8) 


(9) 


(10) 


(11) 


M* Stacks 


7.0 


7.1 




1 














7.1 


7.2 




4 














7.2 


7.3 




4 


14000 ± 600 












7.3 


7.4 




4 


17500 ± 200 












7.4 


7.5 




2 


15700 ± 200 


12800 ± 800 






7.82 ± 0.03 




A^-SFR Stacks 


7.0 


7.1 


0.0 0.5 


1 














7.1 


7.2 


-0.5 0.0 


1 














7.1 


7.2 


0.0 0.5 


2 














7.2 


7.3 


-1.0 -0.5 


1 


13400 ± 500 


11800 ±700 






8.04 ± 0.04 




7.2 


7.3 


-0.5 0.0 


2 















Note. — Column (1): Lower stellar mass limit of the stack. Column (2): Upper stellar mass limit of 
the stack. Column (3): Lower SFR limit of the stack. Column (4): Upper SFR limit of the stack. Column 
(5): Number of galaxies in the stack. Columns (6)-(9): Electron temperatures for [O III], [O II], [N II], 
and [S II]. Column (10): Direct method metallicity. Column (11): N/O abundance. 
(This table is published in its entirety in the electronic edition of the journal. A portion is shown here for 
guidance regarding its form and content.) 



• T e [S II]: -4100 K and -3300 K. 

The T e [O II] and T e [N II] offsets from the T z -T 3 relation for the 
M* stacks are consistent given the scatter, which suggests that 
the T2—T3 relation overestimates the low ionization zone T e 
by ^1000-2000 K. The T e [S II] measurements show a larger 
offset from the T 2 -T 3 relation than T e [0 II] and 7/JN II]. The 
outlier in the M^-SFR 7/ e [S ll]-7/ e [0 III] panel also has a high 
r e [0 II], but this outlier just corresponds to a single galaxy, 
so it may not be representative of all galaxies with this stellar 
mass and SFR. 

The offset between the electron temperatures of the stacks 
and the T2—T3 rela tion is analogous to th e trend for individual 
galaxies found by Pilvugi net alj (J2010), which persists when 
these galaxies are stacked (see Section|4]and Figure^). The 
similar distributions of stacks and individual galaxies relative 
to the 72-73 relation shows that the offset for the stacks is not a 
by-product of stacking but rather a reflection of the properties 
of the individual galaxies (for further discussion see Section 

El. 

At high SFRs (SFR] ;, 5 , and SFRf° 5 ), the offset in r e [0 II] dis- 
appears, and the median T e [0 II] of these stacks is consistent 
with the T2-T3 relation, albeit with a large dispersion. The 
emission from these galaxies is likely dominated by young 
stellar populations, whose hard ionizin g spectrum may be 
similar to the single stellar spectra used by Stasiriska| d 19821) to 
model H II regions. However, a single stellar effective temper- 
ature may not be appropriate for galaxy spectra that include a 
substantial flux contributi on from older H II regions that have 
softer ionizing spectra dKennicutt et al.l l2000b iPilvugin et alj 
l2Ol0l) . 

Figure |4] compares the electron temperatures of [O II] and 
[N II] for the M*-SFR stacks (color-coded by SFR). In the 
two-zone model, both T e [0 II] and r e [N II] represent the 
temperature of the low ionization zone, so these tempera- 
tures should be the same. The stacks scatter around the line 



of equality (black line), though the median offset from the 
r e [0 II] = r e [N II] relation is 1 100 K towards higher r e [N II]. 
If only the stacks that also have detectable [O III] A4363 are 
considered (most of which have r e [0 II] > 8000 K), then the 
median offset from the relation is smaller than the median un- 
certainty on T e [N 11]. The agreement between T e [0 II] and 
r e [N II] for this subset of stacks is consistent with the similar 
offsets found for r e [0 II] and r e [N II] relative to the T 2 -T 3 
relation in Figure [3] 

3.2. Ionic and Total Abundances 

We calculated the ionic abun dance of Q + and O" 
with the nebular.ion ic routine dDe Robertis et al.l 119871: 
IShaw & Dufouri 119951) in IRAF/STSDAS, which determines 
the ionic abundance from the electron temperature, electron 
density, and the flux ratio of the strong emission line(s) rel- 
ative to H/3. We derived the ionic abundance uncertainties 
with the same Monte Carlo simulations used to compute the 
elec tron temperature and density uncertainties (see Section 
13.11 ): the ionic abundance uncertainties were propagated an- 
alytically to calculate the total abundance uncertainties. We 
do not attempt to correct for systematic uncertainties in the 
absolute abundance scale. 

The top two panels of Figure show the ionic abundance 
of + and ++ as a function of stellar mass for the M + stacks 
(open circles) and the M*-SFR stacks (circles color-coded by 
SFR). The + abundance increases with stellar mass at fixed 
SFR and decreases with SFR at fixed stellar mass. The abun- 
dance of ++ is relatively constant as a function of stellar mass 
but is detected in galaxies with progressively higher SFRs as 
stellar mass increases. 

In Figure [5}:, we plot the logarithmic ratio of the O" 1-1 " and 
+ abundances as a function of stellar mass. The dotted line 
in Figure!?}: shows equal abundances of + and O" 1 " 1 ". The con- 
tribution of + to the total oxygen abundance increases with 
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Figure 3. Electron temperatures derived from the [O II], [N II], and [S II] line ratios plotted as a function of electron temperature derived from the [O III] 
line ratio for the M„ stacks (left column) and M*-SFR stacks (right column; color-coded by SFR). The lines in the top, middle, and bottom rows show the 
r c [0 II]-r c [0 III], r c [N II]-r c [0 III], and r c [S II]-r e [0 III] relations (Equation|3), respectively. The outlier in the lower right panel is a single galaxy, so it may 
not be representative of all galaxies with this stellar mass and star formation rate. 



stellar mass at fixed SFR and decreases with SFR at fixed stel- 
lar mass (i.e., in the same sense as how the + abundance 
changes with and SFR). The + abundance dominates 
the total oxygen abundance in the majority of the stacks (i.e., 
above log[M*] = 8.2 for the M* stacks and in half of the M*- 
SFR stacks with detected [O II] AA7320, 7330). Furthermore, 
the + abundance can be measured in many high stellar mass 
and/or low SFR stacks that lack a measured O" 1-1 " abundance, 
which indicates that + is very likely the main ionic species 
of oxygen in these stacks too. A simple extrapolation of the 
log(0 ++ /0 + ) ratio to higher stellar masses for the M* stacks 
shows that the O" 1 ^ abundance would contribute less than 10% 
of the total oxygen abundance. 
We assume that the total oxygen abundance is the sum of 



the ionic abundances of the two dominant species, 

H H+ + H+ ' W 

and the total abundance uncertainties were determined by 
propagating the ionic abundance uncertainties. In highly ion- 
ized gas, oxygen may be found as 3+ , but its contribution to 
the total oxygen abundance is minimal. Abundance studies 
that use the direct method typically measure T e [0 III] and the 
abundance but adopt the T2-T3 relation to infer T e [0 II] 
and the + abundance. However, Figure [3] shows that the T2- 
?3 relation overestimates T e [0 II], which leads to an underes- 
timate of the + abundance and the total oxygen abundance. 
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Figure 4. Electron temperatures derived from the [O II] line ratios as a func- 
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SFR stacks (color-coded by SFR). The line indicates T e [0 II] = 7 C [N II] (as 
assumed in Equation[5J- 
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Figure 5. The ionic abundance of + (panel a), the ionic abundance of 0++ 
(panel b), and the relative ionic abundance of 0++ and + (panel c) as a 
function of stellar mass for the stacks (open circles) and M*-SFR stacks 
(circles color-coded by SFR). The dashed line in panel (c) indicates equal 
abundances of ++ and + . 
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mined from r e [0 II] only (r c [0 III] was inferred with the 7j— T$ relation 
given in Equation|5J and the direct method metallicity determined from both 
r c [0 II] and r e [0 III]. The dotted line denotes the median difference, and the 
dashed line marks the upper mass cutoff for which 7i[0 III] can be indepen- 
dently measured in the M+ stacks. Panel (b): the mass-metallicity relation 
for direct method metallicities determined from T c [0 II] only (gray circles) 
and from both T c [0 II] and T e [0 III] (open circles). To account for the over- 
estimated metallicity (and underestimated T e [0 III]) caused by assuming the 
T2—T3 relation (Equation |3J, we subtract the median metallicity difference 
shown in panel (a) from the T s [0 II]-based metallicities above log(M*) = 9.4 
(shown by the dashed line), which results in the open squares. The arrow 
marks this shift. The sequence of open circles and squares shows the com- 
posite direct method metallicities of the stacks that we will adopt for the 
rest of the paper. We repeated the same procedure for each SFR bin of the 
Mt-SFR stacks. The median metallicity differences are given in Table|4] 

Many of the stacks have measured + and ++ abundances, 
so the total oxygen abundance can be measured accurately in 
these stacks without using the 72-73 relation. 

To extend our total oxygen abundance measurements to 
higher stellar mass, we form a "composite" metallicity cali- 
bration (see Figure |6) that uses the + and O" 1 " 1 " abundances 
when available and the + abundance plus the O" 1 " 1 " abundance 
inferred with the T2-T3 relation if r e [0 II] is measured but 
not r e [0 III] (in the opposite sense from how it is normally 
applied). The total oxygen abundance of the latter group of 
stacks is dominated by the + abundance, so the inferred O" 1 " 1 " 
abundance makes only a small contribution (<10% based on 
the trend indicated by Figure [5];). A simple combination of 
these two metallicity calibrations would lead to a discontinu- 
ity at their interface (in the MZR) because applying the T%— 
73 relation underestimates T e [0 III] and thus overestimates 
the O" 1 " 1 " abundance. To account for this effect, we decrease 
the total oxygen abundances that adopt the T2-T3 relation by 
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the median offset between the two calibrations where they are 
both measured (0.18 dex for the stacks). For the M^-SFR 
stacks, we calculate the median offset for each SFR bin (re- 
ported in Table |4). The offsets are nearly constant as a func- 
tion of and stem from the approximately constant offset 
in the r e [0 ll]-r e [0 III] plot (top row of Figure [3}. Because 
we account for the systematic offset from the 72-73 relation, 
our composite metallicities are insensitive to the exact choice 
of the 72-73 relation. The metallicities of the stacks are pre- 
sented in Table |3] 

Most direct method metallicity studies measure the 
[O III] A4363 line flux but not the [O II] AA7320, 7330 line 
fluxes, so they must adopt a r e [0 ilJ-rjO III] relation, such 
as the 72-73 relation, to estimate the + abundance. One rea- 
son for this is the large wavelength separation between the 
[O II] A3727 strong line and the [O II] AA7320, 7330 auroral 
lines used to measure T e [0 II]. The flux ratio of these two line 
complexes can be affected by a poor reddening correction, 
particularly for the [O II] A3727 line, and some spectrographs 
cannot observe this entire wavelength range efficiently. In in- 
dividual spectra, the [O II] AA7320, 7330 lines can be over- 
whelmed by the noise, which can lead to a lar ge scatter in the 



?a[0 Il]~r e [0 III] diagram (see Fig ure 1 of iKennicutt et alj 
2003 or Figure 4 of Izotov et al. 2006). Fortunately, the noise 
near [O II] AA7320, 7330 appears to be random and is effec- 
tively reduced by stacking, even without the stellar continuum 
subtraction (see Figure|2f). IKennicutt et al. (2003) noted that 
the [O II] AA7320, 7330 line fluxes may be affected by recom- 
bination of O 1 " 1 ", although they find that the typical contribu- 
tion to the [O II] AA7320, 7 330 line fluxes i s <5% (based on 
the correction formulae from lLiu et al.120 00) and that T e [0 II] 
is affected by ~2-3%, which corresponds to <400 K for the 
H II regions in their study. 

We also calculated the ionic abundance of N + with nebu- 
lar.ionic, similar to the procedure used to calculate the ionic 
abundances of + and O^, except that we adopt T e [0 II] 
as the electron temperature instead of T e [N II] because the 
[O II] AA7320, 7330 lines are detected in more stacks and 
with higher SNR than the [N II] A5755 line (see Figure |2). 
The relative ionic abundance of N + /0 + was derived from 
the ionic abund ances of each species. We then assum e that 
N/O = N + /0 + dPeimbert & Costeroll969t lGarnetull990l) to fa- 
cilitate comparison with other studies in the literature (e.g., 
Vila Costas & Edmunds 1993). Although this assumption is 
uncertain, iNava et all (120061) found that it should be accurate 
to -10% for low metallicity objects (12 + log[0/H] < 8.1). 
The N/O abundances of the stacks are reported in Table [3] 

3.3. Strong Line Metallicities 

We compare our direct method metallicities to strong line 
metallicities with various empirical and theoretical calibra- 
tions of the most common line ratios: 

• R 2i : ([O II] A3727 + [O in] AA4959, 5007) / H/3 

• N202: [N n] A6583 / [O n] A3727, 

• N2: [N II] A6583 / Ha, 

• 03N2: ([O in] A5007 / H/3) / ([N n] A6583 / Ha). 

We derived metallicities for our stacks with t he the- 
oretical 7?23 calibrations ofjMcGaugh (119911 here- 
after M91). iZaritsky et al.1 (119941 hereafter Z94), and 
Kobulnick v & Kewle\i (2004 hereafter KK04); the hybrid 
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Figure 7. T e [0 III], r e [0 II], and direct method metallicity for individual 
spectra (small gray cir cles) and stack s in bins of 0.1 dex in stellar mass (large 
black circles) for the Pilyugin et al. (2010) sample relative to the mean of 
galaxies within a stellar mass bin of width 0. 1 dex Mq (shown by the dashed 
line in each panel). The stacks are consistent with the mean r c [0 III], r c [0 II], 
and metallicity within the measurement uncertainties. 



empiri cal-theoretical N2 calibration of iDenicolo et al.l 
(120021 hereafter D02); t he theoretical N202 calibration of 
iKewlev & Dopital ©02, hereafter K D02); and the mostly 
empirical N2 and 03N2 calibrations of lPettini & Pagell (12004 
hereafter PP04). We determined uncertainties on the strong 
line me tallic ities with the Monte Carlo simulations detailed in 
Section lXTl these uncertainties do not account for systematic 
uncertainties in the absolute abundance scale. For a detailed 
discussion of these calibrations and formulae to convert 
between the meta l licitie s derived from each calibration see 
IKewlev & Ellisonl d2008l) . 

4. HOW DOES STACKING AFFECT MEASURED ELECTRON 
TEMPERATURES AND METALLICITIES? 

Stacking greatly increases SNR and thus enables measure- 
ments of physical properties that are unattainable for indi- 
vidual objects. However, measurements from stacked spec- 
tra are only meaningful if they represent the typical prop- 
erties of the objects that went into the stack. To eval- 
uate the effect of stacking on the electron temperatures 
and metallicities, we stacked a sample of 181 SDSS DR6 
dAdelman-McCarfhv et alj|2008l) galaxies with individual de- 
tections of \Q m] A4363 and [O n] AA7320, 7330 from 
Pilvugin et al. (2010) in bins of 0. 1 dex in stellar mass. Figure 
[7] shows the T e [0 III], T e [0 II], and the direct method metal- 
licities of the individual galaxies (gray squares) and stacks 
(black circles) relative to the mean of the galaxies that went 
into each stack. For all three properties, the stacks are con- 
sistent with the mean of the galaxies within the measurement 
uncertainties, which demonstrates that the properties derived 
from galaxies stacked in narrow bins of stellar mass are rep- 
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Figure 8. The electron temperatures derived from the [O II], [N II], and 
[S II] line ratios a s a function of electro n temperature derived from the [O III] 
line ratio for the Pilvuein et al. (2010) sample of galaxies with detectable 
[O III] A4363 and [O II] AA7320, 7330 (squares color-coded by SFR; see 
Section |3J and stacks of the same galaxies in bins of 0.1 dex in stellar mass 
(black circles). The black line shows the T2-T1 relation (Equation[5J. 



resentative of the mean properties of the input galaxies. 

In Figure[8] the [O II], [N II], and [S II] electron tempera- 
tures are plotted as a function of the [O III] electron tempera- 
ture for the galaxies (squares color-coded by SFR) and stacks 
(black circles). The black line in each panel indicates the T2- 
T3 relation (Equation [3). The stacks fall within the distribu- 
tion of galaxies in the T e [0 ll]-7;[0 III] and r e [S ll]-7;[0 III] 
plots (Figure |8k,c). There is some discrepancy between the 
stacks and galaxies in the T e [N ll]-r e [0 III] plot (Figure [8J)), 
but the paucity of [N II] A5755 detections limits the useful- 
ness of any strong conclusions based on T e [N II]. Overall, 
the qualitative agreement between the electron temperatures 
of the stacks and galaxies, especially for r e [0 II] and T e [S II], 
demonstrates that the offset from the T2-T?, relation for the 
stacks shown in Figure [3] is not an artifact of stacking. 

The majority of the galaxies lie below the T2—T3 relation, as 
was previously shown Pilvugi rTet alj (120101). We find a sim- 



ilar result for the galaxies in the T e [S ll]-r e [0 III] relation. 
Galaxies with moderate SFRs (log [SFR] ~ 0.0) are preferen- 
tially further below the T2-T3 relation than galaxies with high 
SFRs (log[SFR] > 1.0) in the r e [0 ll]-7;[0 m l plot. A simi- 
lar effe ct is also present in the M^-SFR stacks. Pilvug irTet alj 
(2010) found that galaxies with lower excitation parameters 
and [O III] A5007/H/3 flux ratios had larger offsets from the 
T2-T3 relation, which is consistent with our result based on 
SFR. They showed that the offset from the T2-T3 relation is 
likely due to the combined emission from multiple ionizing 
sources by comparing the observed r e [0 ll]-r e [0 III] rela- 
tion with the temperature predicted by H II region models 
that include ionizing sources of various temperatures. Based 
on these models, they concluded that differences in the hard- 
ness of the ionizing radiation, caused by the age-dependence 
of H II region spectral energy distributions, govern the scat- 
ter in the T e [0 ll]-r e [0 III] plot for their sample of galaxies. 
Both our results and theirs suggest that galaxies with higher 
SFRs are more similar to the H II region models that served 
as the basis for the T2-T3 relation than galaxies with moderate 
SFRs. This is because they are more likely to be dominated 
by younger stell ar populat i ons that are better approxim ated by 
the input to the Stasihska (1982) models (see Section |3~T1 for 
additional discussion). 

The electron temperatures and metallicities of the stacks are 
unbiased relative to those of the input galaxies, but there is 
some evidence that the integrated galaxy electron temperature 
and metallicity are systematically higher and lower, respec- 
tively, than the electron temperature s and metallicities of the 
individual H II regions in the galaxy. Kobulni cky et alJ (119991) 
compared the electron temperatures and metallicities of indi- 
vidual H II regions in a galaxy to the pseudo-global values 
derived by stacking the spectra of the individual H II regions. 
They showed that the electron temperatures and direct method 
metallicities of their galaxies were biased towards higher tem- 
peratures and lower metallicities by ~ 1000-3000 K and 0.05- 
0.2 dex, respectively, relative to the median values of the 
individual H II regions. Global spectra are biased because 
they are the luminosity-weighted average of the H II regions, 
whose properties can vary widely (see, e.g., the large scatter 
around the 7/2-73 relation f or T e measurements o f individual 
H II regions in Figur e 1 of iKennicutt et al.ll2003l or Figure 4 
of llzotov et al.l 1200a) . The fluxes of the auroral lines might 
be particularly affected by a luminosity-weighted average be- 
cause auro ral fine flux decre a ses no n-linearly with metallic- 
ity. While Ko bulnicky et alj d 19991) only studied the effects 
on [O III] A4363, the relative contribution of each H II region 
likely varies among the commonly measured ionic species, 
potentially yielding results that do not agree with the T2-T3 
relation. We also note that their method of stacking H II re- 
gions does not perfectly simulate global line flux measure- 
ments because it does not account for the contribution of dif- 
fuse ionized gas (i.e., the emission from gas not in H II re- 
gions), which may affect the [N II] and [S II] line fluxes (see 
iMoustakas & Kennicuttl2006T) . In summary, the differences in 
electron temperatures and metallicities between galaxy stacks 
and individual H II regions are dominated by the systematic 
offset between global galaxy properties and individual H II re- 
gions rather than any effects from stacking the global galaxy 
spectra. 

The auroral lines are undetectable in high stellar mass 
galaxies, so we investigate the effect of stacking by com- 
paring the oxygen strong line fluxes of individual galaxies to 
the stack of those galaxies. Figure [9] shows the [O II] A3727 
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Figure 9. [O II] A3727 and [O III] A5007 fluxes relative to H/3 of galaxies 
in one M*-SFR bin (log[M*] = 10.5-10.6 and log[SFR] = 1.0-1.5) and the 
stack of those galaxies. The small black and blue circles represent individual 
galax ies with fluxes measured with our pipeline and the MPA-JHU pipeline 
IT04I) . respectively. The large green circle corresponds to the stack of the 
same galaxies. 



and [O III] A5007 fluxes relative to H/3 for individual galax- 
ies (small black and blue circles) with log(M*) = 10.5-10.6 
and log(SFR) = 1.0-1.5 and the stack of the same galaxies 
(large green circle). The small black and blue circles cor- 
respond to the line fluxes determined with our pipeline and 
the MPA-JHU pipeline, respectively. The distribution of in- 
dividual galaxies with fluxes measured by our pipeline and 
the MPA-JHU pipeline coincide well. In detail, the median 
fluxes from our pipeline are 0.08 and 0.04 dex higher for 
[O III] A5007 and [O II] A3727, respectively, than the median 
fluxes from the MPA-JHU pipeline. The [O III] A5007 and 
[O II] A3727 fluxes of the stack are 0.09 and 0.01 dex higher, 
respectively, than the median of fluxes from our pipeline. Al- 
though the spread is large in the individual galaxies (>1 dex 
for both [O III] A5007 and [O II] A3727), the stack is repre- 
sentative of the typical line fluxes of individual galaxies that 
went into the stack. 

We also note that many of our stacks contain far more 
galaxies than are needed to simply detect a given line, and 
thus are unlikely to be dominated by a few, anomalous galax- 
ies. As an example, we estimate how many galaxies would 
need to be stacked for a detection of [O II] AA7320, 7330. 
If we assume that the uncertainty on the line flux decreases 
as -y/A/gaiaxies, the error on the measurement of any individual 
galaxy is cr stac k * y^Vgdaxies- We use a 5a detection threshold, 
so the minimum number of galaxies needed to detect a line is 
N = [(5er)/flux] 2 . For the M^f-SFRgJ stack, the minimum 
number of galaxies required to detect [O II] AA7320, 7330 is 
^Vgaiaxies = 40, which is well below the actual number of galax- 
ies (1996) in this stack. 

5. THE MASS-METALLICITY RELATION AND 
MASS-METALLICITY-SFR RELATION 

5.1. The Mass-Metallicity Relation 

In Figure [10] we plot the MZR with direct method metal- 
licities for the stacks (circles). We fit the MZR for the M + 
stacks (black line) with the asymptotic logarithmic formula 



suggested by Moustakas e t alJ (120111) : 



12 + log(0/H)= 12 + log(0/H) asm -log f 1+ f \ 1 , 

(5) 

where 12+log(0/H) asm is the asymptotic metallicity, Mjo is 
the turnover mass, and 7 controls the slope of the MZR. This 
functional form is preferable to a polynomial because poly- 
nomial fits can produce unphysical anticorrelations between 
mass and metallicity, particularly when extrapolated beyond 
the mass range over which they were calibrated. The metal- 
licities and fit parameters for the stacks are reported in Tables 
[3]and|4] respectively. For comparison, we show the robust cu- 
bic p olynomial fits of eight str ong line MZRs (col ored lines) 
fro mlKewlev & Ellison! d2008) in Figure [TOk. The IToll IZ94l 
#23 , |KKQ4| fl 2 3 , IKD02J N202, and lM91l flr> M ZRs are based 
on theoretical calibrations, whereas the lD02l N2. |PP04l O3N2. 
and|PP04 N2 MZRs are based on empirical calibrations. In 
Figure [10b . the solid, dashed, and dotted gray lines indicate 
the median, 68% contour, and 95% contour, respectively, of 
the[T03MZR. 

The most prominent aspect of the direct method MZR is 
its extensive dynamic range in both stellar mass and metal- 
licity. It spans three decades in stellar mass and nearly one 
decade in metallicity; this wide range is critical for resolving 
the turnover in metallicity with a single diagnostic that is a 
monotonic relation between line strength and metallicity. The 
broad range in galaxy properties includes the turnover in the 
MZR, which is the first time this feature has been measured 
with metallicities derived from the direct method. Our stacked 
spectra also extend the direct method MZR to sufficiently high 
masses that there is substantial overlap with strong line mea- 
surements, and we use this overlap to compare them. 

The direct method MZR shares some characteristics with 
strong line MZRs but differs in important ways, as can be seen 
in Figure [lOh . The low mass end of the direct method MZR 
starts at log(M^) = 7.4, a full decade lower than the strong line 
MZR s. N onethe less, naive extrapolations of the IT04L ID02L 
IPP04L and[PP04 MZRs are in reasonable agreement with our 
direct method MZR. At a stellar mass of log(M^) = 8.5, the 
lowest stellar mass where strong line MZRs are reporte d, the 
direct method MZR is consistent with the IT04I and the | D02 
MZRs. Above this ma ss, th e direct method MZR and the D02 
MZR diverge from thelT04lMZR. At log(M*) = 8.9, the direct 
method MZR turns over. By contrast, the strong line MZRs 
turns over at a much higher stellar mass (log[MJ ~ 10.5): a 
significant difference that has implications for how the MZR 
is u nders tood in a physical context, which we discuss in Sec- 
tion 17.41 At high mass, the direct method MZR is in good 
agreement with t he empirical stron g li ne cal ibration MZRs, 
but the theoretical lTOll IZ94 IKKOH and lKD02l strong line cal- 
ibration MZRs are offset to higher metallicities by ~0.3 dex at 
log(M + ) = 10.5, the highest mass stack with detected auroral 
lines. Figure [TUb shows the direct method MZR in relation 
to the scatter of the IT04I MZR . The direct method MZR is 
slightly below the median lT04l MZR at log(M*) = 8.5, crosses 
the 16 th percentile at log(M*) = 9.0, and drops below the 2 nd 
percentile at log(M^) = 9.9. Formally, the direct method MZR 
has a scatter of a = 0.05 dex, but this value is not directly 
comparable to the scatter in the IT04I MZR because stacking 
effectively averages over all galaxies in a bin, which erases 
information about galaxy-to-galaxy scatter. 

At low masses (log[MJ = 7.4-8.9; i.e., below the turnover), 
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Figure 10. The direct method mass-metallicity relation for the stacks (circles). In both panels, the thick black solid line shows the asymptotic logarithmic fit 
to the d irect method measurements (see Equation [5). Pane l (a): the colored lines represent various strong line calib rations (Tremonti et al. 2004; Zarits kv et al.l 
ri99llKobulnickv & Kewlev 2004; Kewlev & DopitE 2002 ; McGaugh 1991; Denicolo et al. 2002; Pettini & Page! 2001). Panel (b): the solid, dashed, and dotted 
gray lines show the median, 68% contour, and 95% contour, respectively, of the Tremonti et al. (2004) MZR. The metallicities and fit parameters for the stacks 
are reported in Tables[3]and[4] respectively. 



Table 4 

Mass-Metallicity Relation Fit Parameters 



Stacks 


log(M TO ) 


12+log(0/H) llsm 


7 


Fit Range 


Median log(0/H) Offset 


SFR in Mq yT 1 


[M Q ] 


[dex] 




log(M*) [M Q ] 


[dex] 


(1) 


(2) 


(3) 


(4) 


(5) 


(6) 



MZR 


8.901 


8.798 


0.640 


7.4-10.5 


0.18 


-1.0 < log(SFR) < -0.5 


8.253 


8.726 


0.734 


7.2-9.7 


0.15 


-0.5 < log(SFR) < 0.0 


9.608 


9.118 


0.610 


7.3-10.2 


0.13 


0.0 < log(SFR) < 0.5 


9.836 


8.997 


0.534 


7.6-10.3 


0.12 


0.5 < log(SFR) < 1.0 


27.225 


16.383 


0.449 


8.3-10.6 


0.10 


1.0 < log(SFR) < 1.5 


32.650 


16.988 


0.373 


8.4-10.6 


-0.04 


1.5<log(SFR)<2.0 


28.369 


16.259 


0.438 


9.5-10.8 


0.04 



Note. —Mass-Metallicity Relation given by 12 + log(0/H) = 12 + log(0/H) asra -log(l + (M TO /M*) 7 ). 
Column (1): Stacks included in fits. MZR refers to stacks. The SFR ranges refer to M*-SFR 
stacks. Column (2): Turnover mass. Column (3): Asymptotic metallicity. Column (4): Power-law 
slope. Column (5): Stellar mass range of each fit. Column (6): Median offset between the 
metallicity determined with (i) measured T c [0 II] and infe rred T c [0 III] from the T2-T3 relation and 
(ii) measured T e [0 II] and measured T C [Q III] (see Section [X2t . 



the direct method MZR scales as approximately O/H oc M+ 1 ! 1 . 
While a c ompa rison over the same mass range is not possi- 
ble for the T04 MZR, its low mass slope, as determined from 
log(M„) = 8.5-10.5, is shallower with O/H oc Mj/ 3 . The dis- 
crepan cy in the low mass slopes between the direct method 
and the T04 MZRs could be reasonably explained by the dif- 
ference in the mass ranges over which the slopes were mea- 
sured if the MZ R steepens with decreasing stellar mass (cL, 
iLee et all 12006). We note that the direct method and ID02I 
MZRs have similar slopes and normalizations over a wide 
range in masses from log(M^) = 8.5-10.0. 

5.2. Mass-Metallicity-SFR Relation 

The features of the direct method MZR are shaped by the 
SFR-dependence of the MZR, which we investigate with the 
M + -SFR stacks. Figure QT| shows the M*-SFR stacks (circles 



color-coded by SFR) in the mass-metallicity plane (see Fig- 
ureQ]for the number of galaxies per stack). The solid colored 
lines indicate the asymptotic logarithmic fits (Equation|5]l of 
the M^-SFR stacks of a given SFR, hereafter referred to as 
SFR tracks (e.g., the orange line is the SFRl'j'y track). The 
solid black line is the direct method MZR of the M t stacks 
from Figure [10] the solid, dashed, and dotted gray lines are 
the median, 68% contour, and 95% contour, respectively, of 
the lT04l MZR. The error bars represent the mean error for the 
M*-SFR stacks of a given SFR. 

The M*-SFR stacks help establish the robustness of the di- 
rect method MZR. The low turnover mass and metallicity of 
the direct method MZR relative to the 1T04 and other theo- 
retical strong line calibration MZRs is reminiscent of empir- 
ical strong line calibration MZRs that suffer from a lack of 
sensitivity at high metallicities. However, the most metal-rich 
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Figure 11. The direct method M»-Z-SFR relation for the M*-SFR stacks 
(circles color-coded by SFR) in the mass-metallicity plane. The thick solid 
lines color-coded by SFR show the asymptotic logarithmic fits (see Equation 
|3J for the M*-SFR stacks. The thick black line shows the direct method MZR 
from Figure [lO] The solid, dashed, and dotted gray li nes show the median , 
68% contour, and 95% contour, respectively, of the Tremonti et al. (2004) 
MZR. The error bars correspond to the mean error for the M*-SFR stacks of 
a given SFR. The metallicities and fit parameters for the stacks are given in 
Tables[5]and[4] respectively. 

M^-SFR stacks have some of the highest direct method metal- 
licities (12 + log[0/H] > 9.0) — metallicities well above the 
turnover metallicity of the direct method MZR. These mea- 
surements unambiguously demonstrate that the turnover in the 
direct method MZR is not caused by a lack of sensitivity to 
high metallicities. 

The M*-SFR stacks also can be used to test if galaxies with 
the highest SFRs at a given stellar mass disproportionately in- 
fluence the line fluxes and metallicities of the M± stacks. High 
SFR galaxies have more luminous emission lines and lower 
metallicities and thus may dominate the inferred metallicity 
of the stack. To investigate this possibility, we calculated the 
difference between the metallicity of the M t stack and the 
galaxy number-weighted average of the metallicities of the 
M^-SFR stacks (for the stacks with measured metallicities) at 
a given stellar mass. The median offset is only -0.037 dex 
in metallicity; for reference, the median metallicity uncertain- 
ties for the M* and M*-SFR stacks are 0.019 and 0.027 dex, 
respectively. The slight offset could be due to preferentially 
including the metallicities of M*-SFR stacks with higher SFR 
(lower metallicity) relative to lower SFR (higher metallicity) 
in the weighted average because the former tend to have larger 
line fluxes than the latter, whereas the M* stacks include the 
contribution from galaxies of all SFRs at a given stellar mass. 
Still, the magnitude of this offset is small, which indicates that 
the highest SFR galaxies do not have an appreciable effect on 
the metallicity of the M* stacks because they are quite rare 
(see Figure[TJ. Furthermore, the metallicities of the stacks 
effectively track the metallicity of the most common galaxies 
at a given stellar mass. 

The most striking features of Figure QT|are the 0.3-0.6 dex 
offsets in metallicity at fixed stellar mass between the M*- 
SFR stacks. This trend results from the substantial, nearly 
monotonic dependence of the MZR on SFR. At a given stellar 
mass, higher SFR stacks almost always have lower metallici- 
ties than lower SFR stacks, so there is little overlap between 
the different SFR tracks. Furthermore, the small regions with 



overlap may be the result of the observational uncertainties. 

The interplay between stellar mass, SFR, and metallicity 
for typical galaxies is reflected in the features of the direct 
method MZR, especially the turnover mass. The constant 
SFR tracks (colored lines in Figure Q~T) show that metallicity 
increases with stellar mass at fixed SFR. However, the typical 
SFR also increases with stellar mass, which shifts the "typi- 
cal" galaxy (as measured by the M* stacks) to progressively 
higher SFR and consequently lower metallicity at fixed stellar 
mass. Taken together, the turnover in the MZR is the result 
of the conflict between the trend for more massive galaxies 
to have higher SFRs and the trend for metallicity to decrease 
with SFR at fixed mass. The turnover in the T04 MZR (and 
other strong line calibration MZRs) occurs at a higher stellar 
mass than the the direct method MZR because the strong line 
metallicity calibrations produce a weaker SFR-metallicity an- 
ticorrelation. This means that the progression to higher SFRs 
with increasing stellar mass has less of an effect on the MZR. 

Interestingly, the SFR"„ s stacks (light green circles/line) 
are nearly identical to the T04 MZR in slope, shape, turnover, 
and normalization. While the exact cause of this agreement is 
unclear, it is po ssible that the photoionization models that un- 
derlie the IT04I metallicities assume physical parameters that 
are most appropriate for galaxies with this range of SFR. We 
discuss po tent ial systematic effects of strong line calibrations 
in Section 1731 

The stacks with very high SFRs (SFRj ;g and SFRj j; blue 
and dark blue circles/lines, respectively) have significantly 
lower metallicities than the stack of all galaxies at fixed mass 
in the MZR. The high SFRs and low metallicities of these 
galaxies suggests that they are probably u ndergoing major 
mergers, as found by Peeple s et al.l (120091) for similar out- 
liers. Major mergers drive in considerable amounts of low 
metallicity gas from large radii, which dilutes the metallic- 
ity of the galaxy and triggers vigorous star formation (e.g., 
iKewlev et all 120061 l20iq iTorrev et al.ll20i2h . These stacks 
also have a larger scatter than lower SFR stacks, which is 
likely driven by the small numbers of galaxies per stack cou- 
pled with the large intrinsic dispersion in the individual galaxy 
metallicities. 

5.3. The Fundamental Metallicity Relation 

The orientation of the M^-Z-SFR relation captures 
the im portance of SFR as a se cond parameter to th e 
MZR dLara-Lo pez et al.1 120101: iMannucci etal] 120 1 Oh . 
IMannucci et alj (120101) established the convention that the 
FMR is the projection of least scatter found by choosing a 
free parameter a that minimizes the scatter in the metal- 
licity vs. u n = log(M J - a log(SFR) plane (Equation Q]). 
Mannucci et alJ (120101) found a value of a = 0.32 for a 
sample of SDSS galaxies with me tallicities determine d with 
the semi-empirical calibration of Maio lino et al.l (120081) . As 
metallicity estimates are well known to vary substantially 
between different methods, the parameter a may also be 
different due to potentially different correlations between the 
inferre d metallicity and the SFR. For example, Yat es et al.l 
(120121) used the IT04I metalli cities, rather than those employed 
by IMannucci et alj (120101) . and found a lower value of 
a = 0.19. 

FigureQ~2]shows the fundamental metallicity relation for the 
M*-SFR stacks (circles color-coded by SFR). The scatter in 
metallicity at fixed /i Q is minimized for a = 0.66, which is 
signifi cantl y larger than the a values found by Ma nnucci et alJ 
(2010) and Yate s et al.l (120121) for metallicities estimated with 
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Figure 12. The direct method fundamental metallicity relation for the M*- 
SFR stacks (circles color-coded by SFR). The coefficient (0.66) on log(SFR) 
in the abscissa minimizes the scatter in the FMR (see Equation[T}. The black 
line shows a linear fit to the data, with a slope of 0.43. 



Table 5 

Best Fit a 



Calibration 


a 


(1) 


(2) 


direct method 


0.66 


KK04 


0.24 


M91 


0.17 


Z94 


0.25 


KD02 


0.12 


D02 


0.34 


PP04 N2 


0.30 


PP04 03N2 


0.32 


Note. 


Col- 



umn (1): Metallicity 
cali bratio n (see Sec- 
tion 13.31 for a more 
detailed description 
of the strong line 
calibrations). Col- 
umn (2): The coef- 
ficient on log(SFR) 
in Equation (TJ that 
minimizes the scat- 
ter in the fundamen- 
tal metallicity rela- 
tion. 

strong line calibrations. The scatter for the sta cks differs from 
the sc atte r for individual galaxies (like the iMannucci et alj 
120101 and [Yates et al] 120121 studies) because the number of 
galaxies per stack varies. For a direct comparison, we com- 
puted the value of a for the metallicities derived from vari- 
ous empirical, semi-empirical, and theoretical strong line cal- 
ibrations for the stacks with log(SFR) > -1.0 (the same SFR 
range as the stacks with direct method metallicities) and find 
low a values (a = 0. 12-0 .34) that ar e cons istent with the 
Mann ucci et al) (120101) and lYates et alj (120121) a values (see 
Table |5). The significant difference in a between the direct 
method and the strong line methods indicates that the calibra- 
tions of all of the strong line methods have some dependence 
on physical properties that correlate with SFR. 



The scatter in the direct method FMR (a = 0. 13 dex; Figure 
[12] ) is almost a factor of two smaller than the scatter for the 
M + -SFR stacks with direct method metallicities in the mass- 
metallicity plane (a = 0.22 dex; Figure [TTJ. This decrease is 
due to two features of the M*-SFR stacks at fixed SFR shown 
as the solid colored lines in Figure [TT] (1) they are substan- 
tially offset from one another; (2) they have similar slopes 
with minimal overlap. The former reflects a strong SFR- 
dependence on the MZR; the latter corresponds to a mono- 
tonic SFR-metallicity relation at fixed stellar mass. 

Figure [13] shows the M*-SFR stacks (circles color-coded 
by SFR) in the mass-metallicity plane with metallicities de- 
ter mined with two representativ e strong line calibrations: 
the iKobulnicky & Kewlevl t2004) theo retical /?23 calibration 
(panel a) and the Pettini & Pagell d2004l) empirical N2 calibra- 
tion (panel b). Only stacks with log(M + ) > 8.0 were included 
in Figure[l3]because some stacks at lower stellar masses had 
unphysically high strong line metallicities; to facilitate com- 
parison with Figure QT] only stacks with log(SFR) > -1.0 are 
shown in Figure[T3] The stacks in panel (a) show the metallic- 
ity from the upper branch of j? 23, which were selected to have 
log([N II] A6583/Ha) > -1.1 dKewlev & Elhsodl2008l) . Panel 
(b) shows stacks with -2. 5 < log([N III A6583/H a) < -0.3, 
the calibrated ra nge for th e Pettini & Pagell (12004 N2 calibra- 
tion according to lKewlev"& Ellison (200g). For reference, the 
thick black line shows the direct method MZR. The median, 
68% contour, and 95% contour of the lT04l MZR are indicated 
by the solid, dashed, and dotted gray lines, respectively. 

The scatter in metallicity about the best fit relation de- 
creases only marginally from the MZR to the FMR when 
str ong line cal ibratio ns are used to estimate metallicity. For 
the IKK04I and IPP04I N2 metallicities, the scatter is reduced 
by a = 0.10-^0.09 dex and a = 0.10-S-0.07 dex, respectively. 
Figure [13] also shows that the constant SFR tracks for the 
strong line calibrations in the mass-metallicity plane are both 
more closely packed and overlap more than those of the di- 
rect method. Figure [13] only shows the M^-SFR stacks with 
metallicities from two strong line calibrations, one theoretical 
and one empirical, but the minor reduction in scatter, small 
spread, and considerable overlap are generic features of strong 
line metallicities (the normalization is not). 

A qualitative measure of the spread is the difference be- 
tween the metallicity of the SFR^ 5 (light green) and the 
SFRqj (light blue) stacks at a given stellar mass. There are 17 
stellar mass bins with direct method metallicities for stacks 
with these SFRs. The median metallicity difference for these 
pairs of stacks was 0.38 dex for the direct method, 0.15 dex 
for the Ko bulnickv & Kewlevl d2004l) calibration, and 0.13 dex 
for the lPettini & Pagell d2004l) calibration. The factor of -2-3 
difference between the direct method and strong line metal- 
licities translates into an analogous difference in the SFR- 
dependence of the MZR. 

Another feature of the strong line MZRs at fixed SFR is 
that different SFR tracks turn over at different stellar masses. 
Low SFR tracks turn over at lower stellar masses than high 
SFR tracks, so the sign of the dependence of the MZR on 
SFR changes with stellar mass. At low stellar masses, higher 
SFR stacks have lower metallicities; at high stellar masses, 
the op posite is true — hig her SFR stacks have higher metallic- 
ities. lYates et al] d2012l) found a similar, b ut m ore dramatic, 
result for their sample of galaxies that used |T04| metallicities. 
The origin of the weak SFR-dependence and non-monotonic 
relation for the strong line calibrations is not obvious, but we 
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Figure 13. The M*-SFR stacks (circles color-coded by SFR) in the mass-metallicity plane with metallicities determined with the Kobulnicky & Kewley 1 2004) 
i?23 calibration (panel a) and the Pettini & Pagel 1 2004) N2 calibration (panel b). The thick black l ine shows the direct m ethod MZR from Figure 1 101 The solid, 
dashed, and dotted gray lines show the median, 68% contour, and 95% contour, respectively, of the Tremonti et al. (2004) MZR. The error bars correspond to the 
mean error for the M*-SFR stacks of a given SFR. 
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Figure 14. N + /0 + ratio as a function of direct method oxygen abundance 
(panel a) and M* (panel b) for the M* stacks (open circles) and M*-SFR 
stacks (circles color-coded by SFR). The horizontal lines show the median 
of the low oxygen abundance (12 + log(0/H) < 8.5) and low stellar mass 
(log[M*] < 8.9) data. The positively sloped lines in panels (a) and (b) are 
linear fits to the stacks with 12 + log(CVH) > 8.5 and log(M») > 8.9, respec- 
tively, whose fit parameters are given in Table [6] The error bars show the 
mean error for the stacks (black) and each SFR bin of the M*-SFR stacks 
(color-coded by SFR). If N + /0 + is assumed to trace N/O, as is often done, 
then our results can be compared directly to literature results on N/O. The 
N/O abundances of the stacks are reported in Table [3] 



discuss several potentially relevant effects in Section 1731 

6. N/O ABUNDANCE 

Nitrogen provides interesting constraints on chemical evo- 
lution because it is both a primary and secondary nucleosyn- 
thetic product. The yields of primary elements are inde- 
pendent of the initial metal content of a star but the yields 
of secondary elements are not. In a low metallicity star, 
the majority of the seed carbon and oxygen nuclei that will 
form nitrogen during the CNO cycle are created during he- 



lium burning in the star, so the nitrogen yield of such a 
star will scale roughly with the carbon and oxygen yields. 
In this case, carbon, nitrogen, and oxygen all behave like 
primary elements. After the ISM becomes sufficiently en- 
riched, the nitrogen yield of a star principally depends on 
the amount of carbon and oxygen incorporated in the star 
at birth. The carbon and oxygen still behave like primary 
elements, but nitrogen is a s econdary nucleosynthetic prod- 
uct. Observational studies (IVila Costas & Edmunds 119931; 
van Zee & Havnesl 120061; iBerg et al.ll2012l)~ have found clear 
evidence for primary and secondary n i trogen at low and high 
metallicity. IVila Costas & Edmunds (1993) created a sim- 
ple, closed box chemical evolution model that quantified the 
regimes where nitrogen is expected to behave like a primary 
and secondary element. However, modeling nitrogen enrich- 
ment is difficult because of the large uncertainties in stellar 
yields and the delay time for nitrogen enrichment relative to 
oxygen. Galactic winds also complicate nitrogen enrichment 
becaus e they may preferentially eject oxygen relative to ni- 
trogen (Ivan Zee & H avnes 2006). This is because oxygen is 
formed quickly in massive stars and is available to be ejected 
from galaxies by winds associated with intense bursts of star 
formation. By contrast, the >100 Myr delay before the re- 
lease of nitrogen from intermediate mass AGB stars might be 
sufficient to protect it from ejections by galactic winds. 

In principle, the N/O abundance as a function of oxygen 
abundance can be used to disentangle the effects of nucle- 
osynthesis, galactic inflows and outflows, and different star 
formation histories on the relative enrichment of nitrogen. 
The total N/O ratio is a difficult quantity to measure because 
[N III] lines are not readily observable, so N + /0 + is used fre- 
quently as a proxy for N/O. T his assumption is supported by 
the photoionization models of Garnett ( 1990), which showed 
that the ionization correction factor from N + /0 + to N/O 
should be ^1 to within 20%. Because the ionization factor 
should be close to unity, mos t papers in the literature (e.g., 
IVila Costas & Edmundsl lT993) that show N/O have assumed 
N/O = N + /0 + . For transparency, we plot N + /0 + as a function 



17 



Table 6 

N/O vs. O/H and M* Fit Parameters for M* Stacks 



Abscissa 


Slope 


y-intercept 


Dispersion 


Fit Range 


(1) 


(2) 


(3) 


(4) 


(5) 


12 + log(0/H) 





-1.43 


0.04 


12 + log(0/H) < 8.5 


12 + log(0/H) 


1.73 


-16.15 


0.08 


12 + log(0/H) > 8.5 


M* 





-1.43 


0.04 


M* < 8.9 


M, 


0.30 


-4.04 


0.01 


M* > 8.9 



Note. — Column (1): N/O as a function of 12 + log(0/H) or 
M*. Column (2): Slope of linear fit (set to for first and third rows). 
Column (3): y-intercept of linear fit. Column (4): Dispersion around 
fit. Column (5): Range in 12 + log(0/H) orM* of the fit. 

of direct method oxygen abundance in Figure [14k for the M* 
and M*-SFR stacks. We measured the ionic abundances of 
N + and + with the direct metho d und er the assumption that 
r e [0 II] represents Tz (see Section [3~2l i. 

At low metallicity (12 + log(0/H) < 8.5), we find that the 
M* stacks have an approximately constant value of N + /0 + , 
which is expected for primary nitrogen. These stacks have 
a median of log(N + /0 + ) = -1.43 (indicated by the horizon- 
tal lin^B), which is consi stent with other studies of H I I re- 
gions and dwarf galaxies ( Vil a Costas & Edmunds 1993J)- At 
12 + log(0/H) = 8.5, there is a sharp transition where N + /0 + 
increases steeply with oxygen abundance (slope = 1.73), 
which shows that nitrogen is actin g like a secondary ele- 
ment. Previous observations (e.g., IVila Costas & Edmundsl 
119931) have found a smoother transition between primary and 
secondary nitrogen and a shallower slope in the secondary ni- 
trogen regime, albeit with large dispersion that could be ob- 
scuring these features. The fit parameters of the N + /0 + -0/H 
relation for the M* stacks is presented in Table [6] 

The M* stacks form a tight sequence with a dispersion of 
only a = 0.08 dex, compared to a more t ypical dispersion o f 
a ~ 0.3 dex for individual objects (e.g., iHenry et aD t2000). 
A plausible explanation for the additional scatter in the N/O- 
O/H relation for individual galaxies is the time-dependence 
of N/O caused by the difference in enrichment timescales of 
oxygen and nitrogen following a burst of star formation. The 
M^-SFR stacks show a larger dispersion than the M t stacks, 
potentially because these stacks contain fewer galaxies. The 



SFR%° 5 , and SFR}' 



low and moderate SFR stacks (SFR_j q, 
follow the general trend of the M 4 stacks; however, the high 
SFR stacks (SFR';°, SFRjf,, and SFR?;§) have higher N + /0 + 
at a given oxygen abundance, which may be because these 
galaxies have experienced a large inflow of gas that would 
lower O/H at fixed N/O (i.e., move galaxies to the left in Fig- 
ure fT4h). Another consequence of a vigorous burst of star 
formation is the production of Wolf-Rayet stars that can en- 
rich the gas in nitrogen for a brief peri od before the o xygen 
enrichment from the subsequent SNe II (Ber g et al.l201 lb . We 
see evidence for Wolf-Rayet features, such as He II A4686, in 
some of our stacks, especially at low mass. 

Some of the features in the N + /0 + -0/H relation are clarified 
by the associated relation between N + /0 + and stellar mass, 
which is shown in Figure [14b for the M± and M*-SFR stacks 
(see Table [6] for the fit parameters of the N + /0 + -M* relation 
for the M t stacks). Similar to Figure [14k. there is a primary 

4 We do not show a fit to these points because of the strong a priori ex- 
pectation of a constant N + /0 + at low metallicity (and low mass); however, a 
linear fit would have a slope of —0.21. The analogous slope for the low mass 
N+/0+-M, relation is -0.08. 
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Figure 15. Our direct method MZR (open circles and thick black line) and 
the lLian g et al. 12007) direct method MZR (blue crosses and line). For refer- 
ence, the solid, dashed, and dotted gr ay lines show the me dian, 68% contour, 
and 95% contour, respectively, of the Tremonti et al. (2004) MZR. 

nitrogen plateau in N + /0 + at low stellar mass (log[MJ < 8.9) 
and a steady increase in N + /0 + due to secondary nitrogen en- 
richment above log(M*) = 8.9 (slope = 0.30). However, in 
the secondary nitrogen regime, the N + /0 + -M* relation has a 
much lower dispersion (er = 0.01 dex) than the N + /0 + -0/H 
relation (cr = 0.08 dex). Some of the decreased dispersion 
is due to the larger dynamic range of stellar mass relative to 
oxygen abundance, but the tightness of the N + /0 + -M+ relation 
suggests that the enrichment of nitrogen relative to oxygen is 
well-behaved on average. The essentially zero intrinsic dis- 
persion in the N + /0 + -M* relation can be used to quantify the 
effect of gas inflow and galactic winds on enrichment if all of 
the scatter in N + /0 + at a given O/H is due to gas flows into and 
out of galaxies. As in the N + /0 + -0/H relation, the low and 
moderate SFR stacks (SFRl^, SFR%° 5 , and SFRg jj) roughly 
coincide with the M„ stacks. The high SFR stacks (SFR, 1 ,;", 
SFR[ q, and SFRj ( \) still tend to be more nitrogen enriched 
than the stacks at a fixed M*, but the discrepancy has de- 
creased. The N + /0 + -M* diagram is less sensitive to dilution 
(traced by SFR) because the high SFR galaxies with low O/H 
are less significant outliers when shown as a function of stellar 
mass. 

The N/O-M^ relation ha s been pr eviously investigated 
bvlPerez-Montero & Continil d2009) and lPerez-Montero et al.1 
(2013), who used strong line methods to estimate N/O. They 
found that N/O increased steadily with stellar mass and did 
not show a plateau at low stellar mass associated with primary 
nitrogen enrichment, in contrast to the dir ect method N/O-M* 
relation. However, iPerez-Montero et aD (120 1 3b showed that 
the strong line N/O-M* relation is nearly independent of SFR, 
which is roughly consistent with our finding that the N/O-M* 
relation has only a mild dependence on SFR, particularly at 
log(M*) > 9.0. 

7. DISCUSSION 

7.1. Comparison to a Previous Analysis That Used Auroral 
Lines from Stacked Spectra 

iLiang et ail ((2007) stacked SDSS spectra and applied the 
direct method to estimate the MZR, although their study dif- 
fers from ours in a number of important respects. First, their 
study is based on DR4 spectroscopy of 23,608 galaxies, which 
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is approximately an order of magnitude fewer than our sam- 
ple. Second, they implemented a minimum [O II] A3727 EW 
criterion to select the input galaxies to their stacks in order to 
increase the SNR of their stacked spectra. Finally, they only 
measured r e [0 II] from the [O II] AA7320, 7330 lines and then 
inferred T s [0 III] (and the Q" ionic abund ance) from the Tj— 
T3 relation provided by llzotov et al.l (120061) . These differences 
are likely re s pon sible for the offset between our MZR and the 
ILiang et all (120071) MZR, the absence of a turnover in their 
MZR, and their greater scatter as shown in Figure [T5l 

The [O II] selection criterion can readily expl ain part of 
the offset between our MZRs. Liang et al. (2007) only se- 
lected galaxies with above average [O II] A3727 EW (at fixed 
mass) for galaxies with log(M <r ) < 10 and required a more 
stringent EW([0 II]) > 30 A for galaxies with log(M*) > 10. 
As a result of this selection, their stacks have systematically 
higher SFRs by approximately 0.15 to 0.2 dex. This in turn 
biases the stack s to lower metallicities because of the M+—Z— 
SFR relation dMannucci et al.ll2010HLara-L6pez et al.ll2010l) . 
The magnitude of this effect (^0.05-0.08 dex) accounts for 
part of the difference between the MZRs. Another effect of 
this selection is that the increase in average SFR increases the 
turnover mass and makes it less distinct (see Figure [Til. 

The turnover mass is also not apparent in their MZR due 
to the greater scatter, which is largely due to their order of 
magnitude smaller sample. The scatter around the linear 
fit from log(M+) = 8.0-10.5 for their data is a = 0.12 dex. 
The scatter around an asymptotic logarithmic fit (Equation 
|5]l is reduced only to a = 0.11 dex. An asymptotic loga- 
rithmic fit has an additional degree of freedom relative to a 
linear fit, so the marginal improvement in a suggests that the 
ILiang et all d2007l) MZR can be sufficiently characterized by 
a linear fit. Over the same mass range, the scatter around the 
asymptotic logarithmic fit of our data (thick black line) is only 
a = 0.03 dex, or a factor of four smaller. The smaller scatter 
in our MZR enables a clear identification of th e turnover. 

The method employed bv lLiang et al.l (12007) to estimate the 
oxygen abundance is also distinct from ours and may explain 
the rest of the discrepanc y in the no r mal ization difference be- 
tween our studies. The ILiang et al.l ( 120071) study relies solely 
on the [O II] AA7320, 7330 auroral lines to measure T e [0 II], 
which is used to infer r e [0 III] and the O" 1 " 1 " abundance by ap- 
plying the 7_2-7j re lation and T e [0 lll]-(0 ++ /H + ) formula from 
llzotov et all d2006l) . They did not detect [O III] A4363 in their 
stacks, which they only binned in stellar mass, because they 
had fewer galaxies per bin. The stellar continuum subtraction 
may also have affected the detection of [O III] A4363 because 
of its proximity to the H7 stellar absorption feature, whereas 
the stellar continuum is comparativel y featureless inthe vicin- 
ity of the fO H] AA7320 7330 lines. Liang etafl (120071) used 
the iBruzual & Charlotl (12003b spectral templates, rather than 
the empirical and higher resolution MILES templates that we 
have adopted (see Section |2~3l l, and this difference may also 
have played an important role. As a consequence of their lack 
of a detection of [O III] A4363, their oxygen abundance esti- 
mate depends on the quality o f the assumption th at the galax- 
ies obey the 72-73 relation of llzotov et al.l d20 06). Our empir- 
ical measurements of T2 and T3 indicate that this assumption 
underestimates T3 and overestimates O^/H, which may partly 
explain why our MZRs are in better agreement at high mass 
where + is the dominant ionization state of oxygen. 

7.2. Temperature and Metallicity Discrepancies 



Temperatures and metallicities of H II regions measured 
with the direct method do not always agree with those mea- 
sured with other techniques. For example, temperatures mea- 
sured with the direct method tend to be systematica lly higher 
than those measured from the Balmer continuum (Peimbert 
119671) . Also, the metallicities determined from optical re- 
combination lines (e.g., C II A4267 and O II A4649) and far- 
IR fine-structure lines (e.g., [O III] 52, 88 pm) tend to be 
0.2-0.3 dex higher than those f r om collisionally excited lines 
(Garcfa-Roias & Esteban 2007; Bresolin 2008; Esteban et al. 
2009). The exact cause of these temperature and abundance 
di screpancies is cu rrently not understood. 

Peimbert ( 1967) proposed that temperature fluctuations and 
gradients in H II regions cause direct method temperatures to 
be systematically overestimated, while direct method metal- 
licities are underestimated. To account for temperature varia- 
tions across a nebula, he introduced the concept of t 2 , the root 
mean square deviation of the temperature from the mean. Es- 
timating t 2 has proven to be difficult, so most direct method 
metallicity studies assume f 2 = 0. However, o ptical recom- 
binatio n lines and far-IR fine-structure lines dGarnett et al.l 
I2004al) are less sensitive to temperature than collisionally ex- 
cited lines, so they could be used to estimate t 2 if the discrep- 
ancy between the metallicity determined from collisionally 
excited lines and optical recombination lines or far-IR fine- 
structure lines is assumed to be caused by temperature fluctua- 
tions. The few s tudies that have measured optical recombina- 
tion lines (e.g., iGarcfa-Roias & Esteban 2007; Esteba net al.l 
2009) find that values of f 2 = 0.03-0.07 are necessary to in- 
crease the direct method metallicities by 0.2-0.3 dex to match 
the optical r ecombination l i ne me tallicities. 

Recently, Nic holls et alj d2012l) suggested that the electron 
energy distribution could be the cause of the temperature and 
metallicity discrepancies. Specifically, they questioned the 
widespread assumption that the electrons are in thermal equi- 
librium and can be described by a Maxwell-Boltzmann distri- 
bution. Instead, they suggested that a there might be an excess 
of high energy electrons and proposed that a ^-distribution is 
a more appropriate description of the electron energy distri- 
bution. The ^-distribution is based on direct measurements of 
solar system plasmas. Assuming a ^-distribution for an H II 
region lowers the derived temperature, increases the inferred 
metallicity, and could potentially resolve the discrepancy be- 
tween the temperatures and metallicities found with optical 
recombination lines and collision ally excited lines. 

Models of H II regions by IStasinskal (120051) indicate 
that metallicities based on the direct method could suffer 
from systematic biases in metal-rich H II regions. She 
finds that measuring metallicity from r e [0 III] and T e [N II] 
tends to dramatically underestimate the true metallicity for 
12 + log(0/H) > 8.6 (see her Figure 1). The situation does 
not improve if metallicities are computed with only T e [N II] 
because the derived metallicity can wildly overestimate or 
underestimate the true metallicity depending on the physi- 
cal conditions and geometry of the H II region. However, 
there a re two key differences between the models of Stasihska 
d2005l) and the measurements made in this study that could 
minimize the bias. First, we measured the temperature of the 
low ionization region from T e [0 II], not r e [N II]. Second, we 
analyzed spectra of galaxy stacks and not individual H II re- 
gions, which may average out the large predicted errors. Our 
stacks generally increase smoothly in metallicity as a function 
of M*, which does not rule out systematic error but minimizes 
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the impac t of the individual, catastrophic errors highlighted in 
iStasiriskal (120051) . 

7.3. Strong Line Calibrations and the SFR- dependence of the 

FMR 

The MZRs and FMRs based on strong line calibrations 
(Figure[T3li have a much weaker dependence on SFR than the 
direct method MZR and FMR (Figures [TOifHl ). Relative to 
the direct method MZR and FMR, the strong line MZR and 
FMR have (1) a smaller spread in the mass-metallicity plane 
(compare Figures [TTI and fi"3l. (2) a smaller re duct ion in scat- 
ter from the MZR to the FMR (see Section l5~3l . and (3) a 
smaller value of a (see Table |5j. This trend is a generic fea- 
ture of strong line calibrations that holds for both empirical 
and theoretical calibrations and for all strong line indicators 
(i?23, N2, N202, and 03N2) that we used. Since a strong 
line calibration is only applicable to the physical conditions 
spanned by the calibration sample or model, it is important to 
understand the physical properties of the calibration sample 
for empirical calibrations and the assumptions behind the H II 
region models that underlie theoretical calibrations. 

Figure[T6]compares excitation parameter (P) and R23 (panel 
a) and [O II] and [O III] fluxes relative to H/3 (panel b) for 
galaxies, stacks of galaxies, and H II regions. The gray con- 
tours (50%, 75%, and 95%) and points indicate SDSS star- 
forming galaxies, whos e line flux measurem ents come from 
the MPA-JHU catalog (iTremonfi et alJ 12004 . after we cor- 
rected their measured values for intrinsic reddening. The 
stacks are shown by the open and colored circles. Extragalac- 
tic H II regions with direct method metallicities are repre- 
sented by the light blue contours and crosses. The dereddened 
line flux es of the H II regions c ome from the literature compi- 
lation by Pilv uginetaTl d2012lPl 

Figure [Tol l shows the excitation parameter P as a function 
of /?23- Excitation increases upwards, but/?23 is double-valued 
with metallicity, so metallicity increases to the left for objects 
on the upper branch (the majority of the galaxies and stacks) 
and increases to the right for objects on the lower branch 

5 T he original data ca n be found in Bresolin et al. (2004), Bresol in" et alJ 
(200B1). IBresolirJ (2 007). Bre solin et all (2009a). Bresolin eTaT] (2009bl) 
ICampbell et alj (19 86), Castellanos et al. (2002), de Blok & van der Hulsj 
0993). Estebanetal. (2009), Fierro et al. (1986), French (119801) . 
I Fricke et alj (200ll). [Garnett et al. ilW\). Garnett et all J2004bf) 
iGonzale z-Delgado et al. (1994), Guseva et alj (2000), IGuseva et al] 
420011). IGuseva et alj ( 2003a). Guseva et al. 1 2003b), Gusev a etafj (20041). 
Gusevaetal. (2009), Guseva et al. (2011), Hagele et al. (2008), Hawlev 
(1978), Hodge & Miller 1 1995), Izotov et al. ( 1994), Izotov et al. 1 1997), 
I lzotov & Thuanl (1998al). Ilzoto v & Thuan ~(1998bl). Ilzotov et alj (79991). 
Izotov et al. (2001), Izotov et al. (2004), Izotov & Thuan 12004), Izotov et al. 
(2009), Izotovetal. (2011), Kehria et al. (2004), Kehrig et al. (2011), 
Kennicutt & Skillman 1 2001), Kennicutt et al. 1 2003), Kinkel & Rosa 1 1994), 
Kniazev et alj (2000), iKobulnickv et^!TTl 997), Kunth & Sargent (T9^3I), 
Kwitter & Alleil (19811) ILee et al j (2003al) ILeejiTalT (2003b ). Lee et alJ 
(20041) . ILeeet alj (20051). ILequeux et alj i 19791). IL opez-Sanc hezTTafl 
(2004f). ILopez-Sanchez et alj (20071). ILop ez-Sanchez & Esteban (20091), 
Lopez-Sanchez et al. 12011). Luridiana et al. (2002), Magrini & Goncalves 
(2009), McCall et al. (1985), Melbourne et al. (2004), Melnick et al. 
1 1992), Miller ( 1996), Noeske et al. (2000), Pagel et al. ( 1980), Pagel et al. 
(1992), Pastoriza et al. (1993), Peimbert et al. (1986), Pena et al. (2007), 
Perez-Montero et al. (2009). Popescu & Hopp 12000), Pustilnik et al. 
(20021) . Pustilnik et al. (2003a), Pustilnik etafl (2003bD. IPustilnik et alj 
(200a). Pustilnik et al. (2006), Rav o et alj (T982D. ISaviane et all (2003T 
Sedwick & Aller (1981), Skillman (1985), Skillman & Kennicutt 1 1993), 
Skillman et al. (2003), Stanghellini et al. (2010), Terlevich etajl (75511) . 
Thuan et alj (19951). IThuanet alj (19991) rTorre s-Peimbert etalf 0983) 
Tullmann et al. (2003), van Zee et al. ( 1997), van Zee et al. 1 1998), van Zee 
( 2000), van Zee & Havnes (2006), van Zee et al. (2006), Vilchez et al. 
(1988), Vilchez & Iglesias-Paramo (2003), Webster & Smith (1983), and 
IZahid & Bresolin ( 20Hll 



(most of the compiled H II regions). Figure [16b displays an- 
other projection of the same data in the space defined by the 
dereddened [O II] A3727 and [O III] A5007 line fluxes rela- 
tive to H/3. The dotted lines show constant R23 values, and the 
dashed lines mark constant [O II] A3727/[0 III] A5007 values. 

In Figure[16] the compiled H II regions predominantly over- 
lap with the high excitation and high R23 tail of the galaxy dis- 
tribution in Figure [T6h and the analogous high [O III] A5007 
tail of the galaxy distribution in Figure [T6b . which corre- 
sponds to low metallicity galaxies. The compiled H II re- 
gions have direct method metallicities and therefore at least 
one detectable auroral line, usually [O III] A4363. Because 
the strength of the auroral lines, especially [O III] A4363, is a 
strong function of metallicity and excitation parameter, these 
H II regions were effectively selected to have low metallicities 
and high excitation parameters. Thus, they are not represen- 
tative of the typical conditions found in the H II regions of 
the galaxy sample. Empirical calibrations, which are based 
on samples of H II regions with direct method metallicities, 
are not well constrained in the high metallicity, low excita- 
tion regime where most galaxies and their con stituent H II re- 
gions lie. Forexample, Moustakas et al. (2010) recommended 
only using the empirical Pifyugin & Thuan (2005) R23 cal- 
ibration for objects with P > 0.4. When empirical calibra- 
tions are applied to large galaxy samples, galaxy metallicities 
are systematically undere stimated, particularly at low excita- 
tion and high metallicity (Mousta kas et alj|20l()). Similarly, 
MZRs based on empirical calibrations may have an artificially 
weak dependence on SFR. 

The typical excitation conditions and ^23 values of the 
stacks are much better matched to the overall galaxy distribu- 
tion than the compiled H II regions with direct method metal- 
licities. The stacks probe to both lower excitation (P k 0.2) 
and higher metallicity {R23 ~ 0.4) than the bulk of the com- 
piled H II regions. The stacks do not continue to low R23 
values (<0.3), a region of parameter space populated by the 
most massive and metal-rich galaxies in our sample. The 
[O II] A3727 and [O III] A5007 line fluxes of these galax- 
ies vary significantly, even at the same stellar mass and SFR. 
While the stacks do not reach the lowest R23 values of the 
galaxies, they still trace the average R23 values of the galaxies 
in each stack. 

Theoretical calibrations are based on stellar popula- 
tion s ynthesis models, like STARBURST99 ( iLeitherer et all 
1999), and photoion i zation models, such as MAPPINGS 
( Suth erland &Dopital 119931; iGroves et al l l2004allbl) and 
CLOUDY dFerlandetal.1119981) . The stellar population syn- 
thesis model generates an ionizing radiation field that is then 
processed through the gas by the photoionization model. The 
parameters in the stellar population synthesis model include 
stellar metallicity, age of the ionizing source, initial mass 
function, and star formation history. In the photoionization 
model, the electron density and the ionization parameter are 
adjustable parameters. Because the model grids can span a 
wide range of parameter space, particularly in metallicity and 
excitation parameter, theoretical calibrations have an advan- 
tage over empirical calibrations at high metallicity and low 
excitation, where empirical calibrations are not strictly appli- 
cable. 

However, metallicities derived with theoretical calibra- 
tions can be significa ntly higher (up to 0.7 dex; see 
iKewlev & Ellisonl200 8) than direct method metallicities. The 
most likely cause of this offset is the breakdown of one or 
more of the assumptions about the physics of H II regions 
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Figure 16. Panel (a) shows the excitation parameter, P = [O III] AA4959, 5007 / ([O II] A3727 + [O III] AA4959, 5007), as a function of R^ = ([O II] A3727 + 
[O III] AA4959, 5007) / H/3. Panel (b) shows log([0 II] A3727 / H/3) versus log([0 III] A5007 / H/3). The gray scale contours (50%, 75%, and 95%) and gray 
points correspond to SDSS star-forming galaxies. The white and colored circles represent the M* and M*-SFR (color-coded by SFR) stacks, respectively. The 
light blue contours (50% and 75%) and light blue crosses show H II regions with direct method metallicities from the Pilvuain et al. 12012) compilation. In panel 
(b), the dashed and dotted lines show lines of constant [O II] A3727 / [O III] A5007 and R23, respectively. The stacks trace the overall galaxy distribution better 
than H II regions, especially at lower excitation parameters. The H II regions tend to have high excitation parameters because the auroral line flux is a strong 
function of metallicity and hence R23 . 



in the stellar population synthesis or photoionization mod- 
els. In the stellar population synthesis models, the ionizing 
source is usually treated as a zero age main sequenc e starburst, 
which is not applicable for older star clusters dBerg et al.l 
2011), and the line fluxes can change appreciably as a clus- 
ter (and the assoc i ated H II region) ages. As elucidated by 
iKewle v & Ellison (2008), there are three main issues with the 
photoionization models. First, they treat the nebular geometry 
as either spherical or plane-parallel, which may not be appro- 
priate for the true geometries of the H II regions. Second, 
the fraction of metals deplete d onto dust grains is poorly c on- 
strained by observations (see Draine 2003; Jenkins 2009|) but 
is a required parameter of the photoionization models. Third, 
they assume that the density distribution of the gas and dust 
as smooth, when it is clumpy. While all these assumptions 
might break down to some degree, it is unknown which as- 
sumption or assumptions causes metallicities based on the- 
oretical strong line calibrations to be offset from the direct 
method metallicities, but it is conceivable that the weak SFR 
dependence of theoretical strong line calibration MZRs is also 
due to these assumptions. 

One of t he most intriguing findin gs of the Ma nnucci et alJ 
(2010) and [Lara-Lopez et al.l (f2010) studies is that high red- 
shift observations are consistent with no redshift evolution 
of the strong line FMR up to z = 2.5 and z = 3.5, respec- 
tively. Given the large discrepancies between the local strong 
line and direct method FMRs, a fair comparison between 
the local direct method FMR and a high redshift strong line 
FMR is not possible. An interesting test would be to check 
if high redshift direct method metallicity measurements are 
con sistent with the local direct method FMR. A few stud- 
ies (iHovos et al.l 12005b iKakazu et al J 120071 : 1 Yuan & Kewleyi 
l2009tlErb et al.ll2010tlBra mmer et aLl l2012l) have reported di- 
rect method metallicities at higher redshifts (z ~ 0.7-2.3), 
but none simultaneously provide the stellar masses and SFRs 



of the galaxies. Since the FMR and its evolution provide 
important constraints on theoretical galaxy evolution mod- 
els and form the basis of empirical galaxy evol ution models 
dZahid et al.ll2012bl:lPeeples & Somervillell2013l) . future stud- 
ies that measure all three of these parameters would be valu- 
able. 

7.4. Physical Processes Governing the MZR and M+-Z-SFR 

Relation 

Understanding the baryon cycling of galaxies relies heav- 
ily on the adopted relations between stellar mass, metallic- 
ity, and SFR. Traditionally, the MZR and M^-Z-SFR relation 
have been measured with strong line methods. In this study, 
we have used the more reliable direct method to measure the 
MZR andM^-Z-SFR relation. The direct method MZR (Fig- 
ure [Tojl spans three orders of magnitude in stellar mass from 
log(M*) = 7.4-10.5 and thus simultaneously extends the MZR 
to lower masses by an order of magnitude compared to strong 
line MZRs (e.g., IT04I) and resolves the high mass turnover. 
The features of the direct method MZR that most strongly 
influence the physical interpretations are its low mass slope 
(O/H oc M* 1 / 2 ), its turnover mass (log[MJ = 8.9), and its nor- 
malization (12 + log(0/H) asm = 8.8). The SFR-dependence of 
the MZR (see Figures [TTlandfTZli also serves as an important 
observational constraint for galaxy evolution models. We find 
that the MZR depends strongly on SFR (a = 0.66; Figure fT2l 
at all stellar masses. 

The MZR and M^-Z-SFR relation are shaped by gas in- 
flows, gas outflows, and star formation. The interplay be- 
tween these three proces ses is complex, so hydrodynamic 
galaxy simulations (e.g. , lBrooks et al.ll2007t iFinlator & Pavel 
l2008tlDaye et alJ201 lbtlDave et alJ201 Id) and analytic mod- 
els (e.g.. iPeeples & Shankari 1201 1[ iDave et al.l 120121) have 
been used to establish a framework to interpret the obser- 
vations in a physical context. Below we briefly discuss the 
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physical implications of o ur r esults within the formali sms of 
IPeeples & Shanka ri d20TTh a nd lMnkitor & Pavel d2008b . 

iPeeples & Shankari (1201 ll) developed an analytic model for 
understanding the importance of outflows in governing the 
MZR based on the assumption that galaxies follow zero scat- 
ter relations between stellar mass, gas fraction, metallicity, 
outflow efficiency, and host halo properties. In their for- 
malism, the primary variable controlling the MZR is the 
metallicity-weighted mass-loading parameter, 

> / Zwind \ / Mvjmd \ 

CwM = ) VW ) ' { ) 

where Z w ; nc j and Zism are the wind and ISM metallicities, re- 
spectively, and M w i n d/M* is the unweighted mass-loading pa- 
rameter. ( w i n d can be expressed in terms of the MZR and the 
stellar mass-gas fraction relation by rearranging their Equa- 
tion (20): 

Cwind = 3y Zism - 1 - ctF sas , (7) 

where y is the nucleosynthetic yield, a is a parameter of order 
unity (see their Equation 11), and f gas = M gas /M* is the gas 
fraction. 

If we adopt the Peeples & Shankari d201 11) formalism and 
their fiducial yield and stellar mass-gas fraction relation, then 
we can solve for the M^-Cwind relation implied by the direct 
method MZR. This direct method M^-Cwind relation starts at 
high C wind (Cwind ~ 15) for low mass galaxies (log[MJ = 7.5). 
Then, Cwind decreases with increasing stellar mass, eventually 
flattening and approaching a constant Cwind (C wind ~ 2) above 
the turnover mass (log[MJ = 8.9). Since the D02jMZR has 
a similar shape and normalization to the direct method MZR 
from log(M*) = 8. 5-10 .5, the direct method M*-C w ind rela- 
ti on resembles the ID02I M+- C w ind relation shown in Figure 6 
of IPeeples & Shankari (1201 ll) . Also, the direct method MZR 
implies a similar behavior for Z W md and Zism as a function of 
stellar mass as the ID02I MZR (see their Figure 9). The ratio 
of Zwind/ZisM inversely correlates with how efficiently winds 
entrain ambient ISM. If we adopt the simple relation between 
metallicity and the unweighted mass-loading parameter from 
iFinlator & Day! d2008l) . Zi SM ~ +M wind /M*), then the 
direct method MZR implies an efficiency of mass ejection 

that scales as M w ind/M* oc M~ 1/2 for log(AQ < 9.0. The 
higher Cwind for low mass galaxies relative to high mass galax- 
ies could be due to more enriched winds (larger Z w i n d/ZisM) 
or more efficient mass ejection by winds (larger M w i nd /M*) 
or both. IPeeples & Shankari (120111) found that the M^-Cwind 
relation follows the general shape of the direct method M*- 
Cwind relation regardless of the input MZR (see their Figure 
6). However, the direct method MZR requires more efficient 
metal ej ection by winds than the oretical strong line calibration 
MZRs dT04t IZ94t |KK04 ; IM9 ll) at all stellar masses because 
of the lower normalization of the direct method MZR. We 
note that the yield is poorly constrained, and a higher adopted 
yield requires more efficient outflows to produce the observed 
MZR. 

In contrast to the Peeples & Shankari (1201 ll) framework that 
assumed a zero scatter MZR (and therefore does not account 
for variati ons in the SFR or gas fraction at a fixed stellar 
mass), the Fin lator & Pavel d2008) model, based on cosmo- 
logical hydrodynamic simulations, treats the MZR as an equi- 
librium condition. In their model, galaxies are perturbed off 
the MZR by stochastic inflows but the star formation trig- 
gered by the inflow of gas and the subsequent metal produc- 



tion returns them to the mean MZR. The rate at which galax- 
ies re-equilibrate following an episode of gas inflow sets the 
scatter in the MZR, which is indirectly traced by the SFR- 
dependence of the M^-Z-SFR relation. 

The observed SFR-dependence of the M*-Z-SFR relation 
differs according to the strong line metallicity cali bration used 
to con struct the M*-Z-SFR relation, as found by Yat es et al.l 
(2012) . Specifically, they used metallicities estimated with 
the IMannucci et all d20ld) method and lT04l method. At low 
stellar masses, metallicity decreases with increasing SFR for 
both M*-Z-SFR relations. But at high stellar masses (log[M*] 
> 10.5), the SFR-dependence of the [El M^-Z-SFR rela- 
tion reverses, so th at metallicity i n crease s with increasing 
SFR; however, the IMannucci etall (120101) M*-Z-SFR rela- 
tion collapses to a sin gle sequence that is independent of 
SFR. I Yates et ail d2012l) su ggested that the SFR-dependence 
of the lMannucci et alJd20loh M^-Z-SFR relation at high stel- 
lar mass is obscured by the N2 indicator (which was aver- 
aged with the metallicity estimated from R23) used in the 
IMannucci et alj d2010l) metallicity calibration, which satu- 
rates at high metallicity. 

Unlike the IMannucci et alj ( 120101 ) and[j03M*-Z-SFR re- 
lations, the SFR-dependence of the direct method M*-Z-SFR 
relation does not change dramatically with stellar mass. There 
is little overlap between the constant SFR tracks in the di- 
rect method M^-Z-SFR relation (Figure [PTV Furth erm ore, 
the SFR-dependence is strong (a = 0.66; see Section IB31 l. so 
the scatter in the direct method MZR for individual galaxies 
(if it could be measured) w oul d be larger than the scatter in 
the IMannucci et al.l (12010) and IT04I MZRs. Within the con- 
text of the Finlator & Dave (2008) model, this means that the 
direct method MZR i mplies a longer timesca le fo r gala xies to 
re-equilibrate than the lMannucci et al.l (120101) and|T04 MZRs. 
We note that the direct method M^-Z-SFR relation does not 
probe above log(M^) = 10.5 because the auroral lines are un- 
detected in this regime; how ever, this mass scale is w here 
the discrepancies between the Man nucci et al.l (120101) and T04 
metallicities are the largest — potentially due to a break down 
of s trong line calibrations at high metallicities (see Section 

E3. 

8. SUMMARY 

We have measured [O III], [O II], [N II], and [S II] electron 
temperatures, direct method gas-phase oxygen abundances, 
and direct method gas-phase nitrogen to oxygen abundance 
ratios from stacked galaxy spectra. We stacked the spectra 
of ~200,000 SDSS star-forming galaxies in bins of (1) 0.1 
dex in stellar mass and (2) 0.1 dex in stellar mass and 0.5 
dex in SFR. The high SNR stacked spectra enabled the de- 
tection of the temperature-sensitive auroral lines that are es- 
sential for metallicity measurements with the direct method. 
Auroral lines are weak, especially in massive, metal-rich ob- 
jects, but we detect [O III] A4363 up to log(M*) = 9.4 and 
[O II] AA7320, 7330 up to log(AQ = 10.5, which is gener- 
ally not feasible for spectra of individual galaxies. We used 
the auroral line fluxes to derive the [O III] and [O II] electron 
temperatures, the ++ and + ionic abundances, and the total 
oxygen abundances of the stacks. 

We constructed the direct method mass-metallicity and 
M*-Z-SFR relations across a wide range of stellar mass 
(log[MJ = 7.4-10.5) and SFR (logfSFR] = -1.0-K2.0). 
The direct method MZR rises steeply (O/H oc M* 1 / 2 ) from 
log(M*) = 7.4-8.9. The direct method MZR turns over at 
log(Mi) = 8.9, in contrast to strong line MZRs that typi- 
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cally turn over at higher masses (log[MJ ~ 10.5). Above 
the turnover, the direct method MZR approaches an asymp- 
totic metallicity of 12 + log(0/H) = 8.8, which is consistent 
with empirical strong line calibration MZRs but ~0.3 dex 
lower than theoretica l strong line calibration MZRs like the 
iTremonti et all (S2004) MZR. Furthermore, we found that the 
SFR-dependence (as measured by the value of a that mini- 
mizes the scatter at fixed \i a = log(M*) - alog(SFR) in the 
fundamental metallicity relation; see Equation[f]i of the direct 
method M*-Z-SFR relation is ^2-3 times larger {a = 0.66) 
than for strong line M^-Z-SFR relations (a ~ 0.12-0.34). Its 
SFR-dependence is monotonic as a function of stellar mass, 
so constant SFR tracks do not overlap, unlike strong line M*- 
Z-SFR relations. 

We also showed that the direct method N/O relative abun- 
dance correlates strongly with oxygen abundance and even 
more strongly with stellar mass. N/O exhibits a clear transi- 
tion from primary to secondary nitrogen enrichment as a func - 
tion of oxygen abundance and stellar mass. 

The slope, turnover, normalization, and SFR-dependence 
of the MZR act as critical constraints on galaxy evolution 
models and are best measured by methods that do not rely 
on strong line diagnostics, such as the direct method. Future 
work should aim to construct a direct method MZR of indi- 
vidual galaxies with high SNR optical spectra that enable the 
detection of auroral lines in high mass and high metallicity 
objects. Furt hermore, metalliciti es based on Herschel Space 
Observatory (Pilbra tt et al.ll2010l) measurements of the far-IR 
fine-structure lines (Croxall et al. in prep.) from the KING- 
FISH survey dKennicutt et alj|201 11) will provi de a valuable 
check on the absolute abundance scale (see also lGarnett et al.l 
2004a), which is a major outstanding uncertainty for galaxy 
evolution studies. These types of investigations will improve 
our understanding of the galaxy formation process, particu- 
larly the cycling of baryons between galaxies and the IGM. 
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